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PREFACE 

The  work  reported  herein  was  conducted  by  the  Arnold  Engineering  Development  Center  (AEDC),  Air 
Force  Materiel  Command  (AFMC),  under  funding  from  the  Air  Force  Office  of  Scientific  Research 
(AFOSR).  The  results  of  the  research  were  obtained  by  Aerospace  Testing  Alliance,  the  operations, 
maintenance,  information  management,  and  support  contractor  for  AEDC,  AFMC,  Arnold  Air  Force 
Base,  Tennessee,  under  Project  Number  9830.  The  Air  Force  Project  Managers  were  Drs.  Neil  Glassman, 
Thomas  Beutner,  and  John  Schmisseur  for  AFOSR  and  Mr.  Steven  A.  Bancroft  for  AEDC.  The  author  is 
also  indebted  to  his  colleagues  at  AEDC,  in  particular,  Mr.  Marvin  Sellers  and  Mr.  Daryl  Sinclair. 
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1.0  INTRODUCTION 

The  Arnold  Engineering  Development  Center  (AEDC)  is  one  of  the  leading  aerospace  test  centers  in 
the  world  and  operates  many  facilities  that  have  unmatched  capabilities  anywhere.  In  addition,  AEDC 
strives  to  offer  its  customers,  both  military  and  commercial,  the  best  measurement  techniques  available. 
Among  the  newest  of  these  techniques  is  the  use  of  pressure-sensitive  paint  (PSP)  for  the  measurement  of 
full-field  pressure  distributions  on  the  surfaces  of  test  articles. 

Introduced  in  the  early  1990s,  PSP  techniques  have  since  gone  through  a  significant  maturation 
process,  but  prior  to  the  start  of  Air  Force  Office  of  Scientific  Research  (AFOSR)  contract  F40600-03-R- 
0001,  a  principal  challenge  remained:  automated,  online  processing  of  the  large  amounts  of  image  data 
that  are  acquired  during  a  typical  PSP  test.  It  was  not  unusual  for  such  data  reduction  to  take  weeks 
following  completion  of  the  test.  This  was  unsatisfactory  to  test  customers  and  has  posed  a  significant 
challenge  to  scientists  and  engineers.  This  challenge  was  the  motivation  for  the  present  AFOSR  contract. 

This  report  shows  that  the  image-processing  challenge  has  been  met  successfully,  in  large  part 
because  of  the  support  of  AFOSR.  This  success  is  rooted  not  merely  in  the  development  of  improved 
software,  per  se.  Several  noteworthy  innovations  have  been  developed  with  regard  to  automatic  image 
registration  and  parallel  processing — innovations  that  have  the  potential  for  leading  to  improved  data- 
processing  techniques  in  other  disciplines. 

Furthermore,  as  a  direct  result  of  the  support  provided  by  AFOSR,  opportunities  have  arisen  to 
improve  not  only  the  speed  with  which  data  from  a  PSP  test  can  be  processed,  but  also  the  quality  of  the 
resulting  data.  In  addition,  several  studies  have  been  performed  to  guide  the  continued  development  of 
improved  paints  and  data  acquisition  strategies  for  the  application  of  lifetime-based  PSP  techniques  (for 
example,  the  kind  of  strategies  being  pursued  under  AFOSR  contract  STTR-AF04-T001,  “Improved 
Pressure-  and  Temperature-Sensitive  Paint.”) 

This  report  consists  of  previously  published  material  that  was  made  possible  with  AFOSR  support 
under  this  contract.  The  material  has  been  reformatted  into  a  cohesive  report,  with  cross-references 
between  sections,  where  applicable.  A  complete  list  of  publications  generated  under  this  contract  follows, 
with  a  brief  explanation  of  how  this  material  was  merged  into  the  present  report. 

SOURCE  DOCUMENTS 

1.  W.  Ruyten  and  M.  Sellers,  “On-Fine  Processing  of  Pressure-Sensitive  Paint  Images,”  Journal  of 
Aerospace  Computing,  Information,  and  Communication,  Vol.  1,  September  2004,  pp.  372-382.  This 
material  is  presented  first  (in  Section  2.0  of  the  present  report),  because  it  provides  a  good 
introduction  both  to  the  problem  of  automated  processing  and  to  the  approach  to  its  solution. 

2.  W.  Ruyten  and  M.  Sellers,  “On-Fine  Processing  of  Pressure- Sensitive  Paint  Images,”  AIAA  Paper 
2003-3947,  21st  AIAA  Applied  Aerodynamics  Conference,  Orlando,  Florida,  June  23-26,  2003. 
This  paper  is  an  earlier  version  of  the  material  in  Doc.  1,  above. 

3.  W.  Ruyten  and  W.  E.  Sisson,  “Message-Passing  for  Parallel  Processing  of  Pressure-Sensitive  Paint 
Images,”  DoD  HPCMP  Users  Group  Conference,  Williamsburg,  Virginia,  June  7-11,  2004, 
published  by  IEEE  as  ISBN  0-7695-2259-9,  pp.  308-312.  This  material  is  presented  in  Section  3.0  of 
the  present  report. 

4.  W.  Ruyten,  “Automatic  Image  Registration  for  Optical  Techniques  in  Aerodynamic  Test  Facilities,” 
AIAA  Paper  2004-2400,  24th  AIAA  Aerodynamic  Measurement  Technology  and  Ground  Testing 
Conference,  Portland,  Oregon,  June  28-July  1,  2004  [one  of  two  papers  (out  of  a  field  118  papers)  to 
be  presented  with  an  outstanding  paper  award].  This  material  is  presented  in  Section  4.0  of  the 
present  report . 
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5.  W.  Ruyten,  “Sting  Deflection  Model  for  Automatic  Image  Registration  in  Wind  Tunnel  Testing,” 
Journal  of  Aircraft,  Vol.  43,  No.  3,  May  -  June  2006,  pp.  853-855.  This  material  was  extracted  from 
Doc.  4,  above. 

6.  W.  Ruyten,  “Implementation  of  Template-Based  Image  Registration  for  Aerodynamic  Testing,” 
Journal  of  Aerospace  Computing,  Information,  and  Communication,  Vol.  3,  No.  4,  April  2006,  pp. 
138-145.  This  material  also  was  extracted  from  Doc.  4,  above. 

7.  W.  Ruyten,  “Extension  of  Subpixel  Augmented  Template  Matching  to  Imperfectly  Rendered  Binary 
Targets,”  Optical  Engineering,  Vol.  43,  No.  3,  March  2004,  pp.  639-647.  This  paper  was  initiated 
prior  to  the  AFOSR  contract  as  an  IEEE  conference  paper  but  was  brought  to  completion  as  a 
refereed  publication  with  the  support  of  AFOSR.  The  material  is  presented  in  Section  5.0  of  this 
report. 

8.  W.  Ruyten  and  M.  Sellers,  “Improved  Data  Processing  for  Pressure- Sensitive  Paint  Measurements  in 
a  Production  Facility,”  AIAA  Paper  2006-1042,  44th  Aerospace  Sciences  Meeting  and  Exhibit, 
Reno,  Nevada,  9-12  January,  2006.  This  material  forms  the  culmination  of  the  work  that  was 
performed  under  the  AFOSR  contract.  It  reflects  the  innovations  that  were  made  possible  with  the 
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2.0  ON-LINE  PROCESSING  OF  PRESSURE-SENSITIVE  PAINT  IMAGES* 

2.1  OVERVIEW 

Section  2.0  documents  a  milestone  that  has  been  reached  at  AEDC,  namely  to  reduce  data  reduction 
times  for  both  intensity-based  and  lifetime-based  pressure-sensitive  paint  measurements  from  an  eight- 
camera  system  sufficiently  that  one  data  point  can  be  processed  in  full  (in  5  to  1 0  seconds)  while  the  next 
data  point  is  being  acquired.  This  paper  describes  the  hardware,  software,  and  data  reduction  strategies 
used,  and  it  presents  results  for  benchmark  data  from  wind  tunnel  tests  on  scale  models  of  an  F-16C 
fighter  jet  and  NASA’s  X-38  Crew  Return  Vehicle. 

2.2  INTRODUCTION 

Pressure-sensitive  paint  (PSP)  has  established  itself  as  an  important  test  and  evaluation  tool  in  the 
area  of  aerodynamic  testing,  particularly  for  use  in  transonic  wind  tunnels. Most  commonly,  the 
luminescence  from  the  painted  test  article  is  captured  using  digital  cameras.  Thus,  PSP  is  essentially  an 
image-based  technique.  (The  much  less  frequently  used  approach  based  on  a  scanned  laser  spot  system  is 
not  considered  here.)  To  achieve  full  coverage  across  the  surface  of  the  test  article,  multiple  cameras 
typically  are  used.  Data  are  then  processed  and  presented  to  the  test  customer  either  in  2D  image  format 
or  as  data  that  have  been  mapped  to  a  3D  grid  of  the  test  article.  The  latter  format  is  preferred  because  it 
allows  force  and  moment  integrations  to  be  performed,  across  either  the  whole  surface  or  parts  thereof 

Considerable  computational  resources  are  required  to  convert  raw  images  into  fully  processed  PSP 
data.  In  the  past,  days  or  weeks  were  required  to  perform  this  task.  Over  the  years,  analysis  techniques 
have  become  more  sophisticated  and  computers  have  become  faster,  making  it  possible,  at  present,  to 
process  data  while  a  wind  tunnel  test  is  in  progress.  The  purpose  of  this  section  is  to  review  the 
technology  that  has  made  this  milestone  possible  at  AEDC. 

2.3  THEORY 

2.3.1  Paint  Response 

The  theory  behind  PSP  has  been  reviewed  extensively  elsewhere.’’^  At  the  heart  of  PSP  processing  is 
the  notion  that  the  pressure,  P,  can  be  obtained  from  a  ratio,  R,  of  two  signals,  and  by  a 
calibration  function,  P{R,  T),  where  T  is  the  temperature.  Symbolically,  this  may  be  written  as 


and 

Pij^P{Rij,T).  (2) 

The  subscripts  ij  indicate  that  the  respective  quantities  may  be  thought  of  as  pixel  values  in  a  digital 
image.  The  pressure  from  Eq.  (2)  is  usually  converted  to  a  pressure  coefficient  according  to  the  relation 

Cp,ij-iPij-Poo)/Q,  (3) 

where  Poo  and  Q  are  the  freestream  static  and  dynamic  pressures,  respectively. 


W.  Ruyten  and  M.  Sellers,  “On-line  Proeessing  of  Pressure-Sensitive  Paint  Images,”  Journal  of  Aerospace 
Computing,  Information,  and  Communication,  Vol.  1,  September  2004,  pp.  372-382.  An  earlier  version  was 
presented  as  AIAA  Paper  2003-3947,  21st  AIAA  Applied  Aerodynamies  Conferenee,  Orlando,  Florida,  June  23-26, 
2003. 
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In  intensity-based  PSP  tests,  the  signals  5^’^  and  5^^'  in  Eq.  (1)  are  the  wind-off  and  wind-on  signals.  In 
lifetime-based  PSP  tests,  both  5^'^  and  are  obtained  at  the  run  condition,  but  at  different  gates  with 
respect  to  the  excitation  pulse  train. 

The  calibration  function  P{R,T)  typically  is  cast  in  the  form  of  a  polynomial  with  fixed  coefficients.^’^ 
Figure  1  shows  an  example  of  a  set  of  calibration  curves  for  the  paint  PtTFPP  in  FIB,^  which  is  produced 
by  Innovative  Scientific  Solutions,  Inc.  and  is  the  current  industry  standard  for  large-scale  PSP  testing  in 
the  U.S.  Ideally,  the  temperature,  T,  in  Eq.  (2)  would  be  provided  on  a  per-pixel  basis.  In  practice,  a 
single,  average  value  typically  is  supplied. 


Figure  1.  Example  of  PSP  Calibration  Curves 

More  complex  implementations  of  PSP  are  possible — for  example,  those  based  on  the  use  of 
bilumiphore  paints  (Refs.  10  and  11)  or  those  based  on  multigate  imaging.'^"'^  In  such  alternate  schemes 
two  signal  ratios  are  defined,  allowing  the  determination  of  either  a  temperature-corrected  PSP 
measurement  or  a  combined  pressure  and  temperature  measurement. 

In  the  following  paragraphs  it  is  assumed  that  the  paint  response  is  given  by  Eq.  (2)  and  that  the 
calibration  function,  P{R,  T),  has  been  determined  prior  to  the  wind  tunnel  test. 

2.3.2  Image  Mapping 

To  establish  the  mapping  transformation  between  2D  image  coordinates  and  3D  spatial  coordinates 
on  the  surface  of  the  test  article,  registration  targets  are  applied  to  the  test  article.  Eet  X„  denote  the  3D 
model  coordinates  of  the  target  with  index  n,  as  determined  by  a  coordinate-measuring  machine.  The 
resulting  image  coordinates,  (x„,y„),  may  be  written  symbolically  as  (see  Section  4.5)''*’'^: 

’ T«  )  ~  ’Pc  > Qc ’ ’  t)>  (4) 

where  the  function  represents  a  projective  transform  that  takes  as  input  the  3D  model  coordinates, 
X„,  of  the  targets;  the  set  of  exterior  parameters,  pc*,  for  camera  c  (representing  the  position  and 
orientation  of  camera  c  with  respect  to  the  test  facility);  the  set  of  interior  parameters,  Qc,  for  the  same 
camera  (representing  the  effective  focal  length,  optical  center,  and  lens  distortion  parameters);  a  set  of 
angles,  a,„,  that  describes  the  pitch  and  roll  settings  of  the  model-positioning  system  at  a  model  attitude  w; 
and  a  fixed  set  of  model  alignment  parameters,  t. 

Before  the  function  F{...)  can  be  determined  for  a  given  image,  the  image  coordinates,  {x„,yr),  of  a 
sufficiently  large  set  of  targets  must  be  obtained.  This  requires  that  the  registration  targets  be  located  and 
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identified  in  the  image  (see  Fig.  2).  To  do  this  automatically  and  efficiently  has  been  the  single  biggest 
challenge  in  achieving  on-line  processing  of  PSP  data.  Aspects  of  this  task  are  described  in  Refs.  16 
through  18  and  in  Section  4.0  of  this  report. 


Figure  2.  Example  of  PSP  Image  with  Registration  Targets 

Once  the  mapping  function,  F,  is  known,  it  may  be  used  to  map  image  pixels  to  the  predefined  3D 
grid  of  the  test  article.  This  typically  involves  the  use  of  a  z-buffer  to  preclude  image  data  from  being 
mapped  to  occluded  parts  of  the  3D  surface  grid.  So  far,  only  multizone  structured  grids  have  been  used 
as  3D  grids,  though  the  use  of  unstructured  grids  is  under  consideration. 

A  related  image  registration  problem  in  intensity-based  PSP  is  that  of  aligning  a  wind-on  image  to  a 
wind-off  image.  If  the  misalignment  between  the  two  images  is  not  too  severe,  the  image  coordinates  in 
the  second  image  may  be  related  to  the  image  coordinates  in  the  first  image  by  a  polynomial 
transformation  that  may  be  written  symbolically  in  terms  of  a  function,  G,  as*®’^**: 

=  (5) 

2.3.3  Self-Illumination  (SI)  Correction 

When  emitted  light  from  a  low-pressure  region  on  the  surface  of  the  test  article  (i.e.,  a  region  with 
relatively  high  signal)  is  emitted  toward  a  high-pressure  region  (i.e.,  a  region  with  relatively  low  signal) 
and  is  reflected  toward  the  camera,  the  increased  signal  caused  by  the  reflection  process  tends  to  lower  the 
measured  pressure.  To  a  lesser  extent,  measured  pressures  in  low-pressure  regions  are  affected  as  well.  It 
is  possible  to  correct  for  these  effects  by  approximating  the  painted  surface  as  a  diffuse  reflector.  The  net 
emitted  signal,  Sk,  at  a  point,  k,  on  the  3D  grid  can  then  be  obtained  from  the  measured  signals.  Si,  at 
other  points,  /,  according  to 

S'k  -  Rdiff  X  ^kl^l  ’ 

I  (6) 

where  Rdiff  is  the  reflectivity  of  the  paint,  is  an  influence  coefficient  that  depends  only  on  the  geometry 
of  the  test  article,  and  the  summation  is  over  that  part  of  the  surface  with  nonzero  influence  coefficients.^' 
Higher-order  corrections,  based  on  the  use  of  a  bidirectional  reflection  distribution  function,  are  also 
possible.^^ 
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The  SI  correction  from  Eq.  (6)  needs  to  be  performed  separately  on  both  of  the  signal  terms  in  Eq. 
(1).  This  implies  that  the  image  data  must  be  mapped  to  the  3D  grid  before  a  signal  ratio  is  calculated.  In 
this  regard,  SI  corrections  require  a  deviation  from  the  conventional  data-processing  sequence  in  which 
data  from  each  camera  are  processed  separately  to  a  pressure  distribution  on  the  3D  grid. 

The  A-matrix  from  Eq.  (6)  can  be  precalculated,  allowing  the  SI  correction  to  be  performed 
efficiently.  Moreover,  because  the  reflected  light  is,  by  assumption,  diffusely  distributed,  the  SI  correction 
needs  to  be  calculated  explicitly  only  for  a  fraction  of  the  grid  points  (10  percent,  typically).  The  SI 
correction  at  the  remaining  points  can  then  be  obtained  by  interpolation. 

2.4  PROCESSING  SCHEMES 

More  than  70  macro  commands  are  built  into  the  original  image-processing  software  (that  is,  the 
Green  Boot  code;  see  Section  2.6.1).  The  most  important  of  these  can  be  grouped  into  the  following  eight 
processing  steps.  (These  eight  reflect  more  closely  the  processing  commands  in  the  more  recently 
developed  software.) 

DEFINE:  Eocate  the  required  input  data  files  on  the  basis  of  the  to-be-processed  run  and  sequence 
numbers,  and  generate  the  names  of  the  resulting  output  files. 

EOAD:  Eoad  data  files  into  memory  and  perform  one  or  more  of  the  following:  Subtract  a  black 
image  (an  “image”  collected  with  the  lens  covered);  divide  by  a  flat-field  image  [an  “image”  collected 
without  a  lens  on  the  camera,  with  the  charge-coupled  device  (CCD)  array  exposed  to  a  uniform  light 
field];  scale  the  image  by  an  exposure  time  or  by  the  number  of  pulses  over  which  the  integration  was 
performed;  and/or  apply  a  threshold  value  below  which  image  data  are  considered  to  be  indistinguishable 
from  the  background.  Because  of  the  use  of  image  ratios  in  Eq.  (1),  flat-field  corrections  are  usually  not 
necessary. 

REGISTER:  Find  (to  subpixel  precision)  and  identify  all  of  the  visible  registration  targets  in  the 
image  and  determine  the  resulting  mapping  transformation  between  3D  model  coordinates  and  2D  image 
coordinates,  as  expressed  by  Eq.  (4).  Two  tasks  typically  performed  along  with  image  registration  are  the 
calculation  of  the  location  of  visible  pressure  taps  (for  use  with  in-situ  calibrations)  and  the  removal  of 
registration  targets  and  other  artifacts  by  patch-interpolation. 

AEIGN:  Given  the  image  coordinates  of  registration  targets  in  a  set  of  wind-off  and  wind-on  images, 
determine  the  coefficients  of  the  image  alignment  transformation  from  Eq.  (5)  and  apply  the 
transformation  to  align  the  two  images.  This  alignment  is  only  needed  for  intensity-based  PSP 
measurements  and  is  to  compensate  for  the  shift  in  model  position  that  occurs  as  a  result  of  sting 
deflection. 

CONVERT:  Calculate  the  ratio  from  Eq.  (1),  either  between  two  images  or  between  two  sets  of  data 
on  the  same  3D  grid.  Set  the  result  to  “undefined”  for  points  that  do  not  have  a  valid  value  in  both  sets; 
then  convert  the  ratios  to  pressure  using  Eq.  (2),  or  to  a  pressure  coefficient  using  Eq.  (3).  If  an  in-situ 
calibration  is  performed,  the  paint-derived  pressure  values  are  offset  so  as  to  achieve  optimum  agreement 
between  the  PSP  values  in  the  immediate  vicinities  of  the  pressure  taps  and  the  tap  values  themselves. 

MAP:  On  the  basis  of  Eq.  (4),  map  2D  image  data  to  the  predefined  3D  grid  of  the  test  article.  To 
prevent  mapping  to  occluded  parts  of  the  3D  grid,  a  z-buffer  algorithm  is  used  to  determine  visibility.  If 
the  same  grid  point  is  mapped  by  two  or  more  cameras,  the  camera  with  the  view  most  nearly  normal  to 
the  surface  is  used  or  an  average  value  is  computed,  as  is  described  in  Section  6.4. 

SI-CORRECT:  Perform  the  self-illumination  correction  from  Eq.  (6). 

SAVE:  Save  intermediate  or  final  data  to  the  hard  disk  in  3D  grid  format.  In  addition  to  the  pressure 
or  pressure  coefficient  at  every  grid  point,  the  saved  file  may  contain  additional  values,  such  as  the 
calculated  uncertainties  of  the  pressure  data  (See  Section  6.5). 
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Figure  3  shows  two  processing  sequences  based  on  the  eight  processing  steps  above.  Both  sequences 
lead  from  two  sets  of  raw  image  data  (sets  “1”  and  “2”  in  Fig.  3)  to  pressure  coefficients  on  a  3D  grid  (set 
“3”  in  Fig.  3).  The  same  steps  are  followed  for  both  intensity-based  and  lifetime -based  PSP.  The  primary 
difference  between  the  two  sequences  is  the  point  at  which  2D  image  data  are  mapped  to  the  3D  grid.  In 
the  second  sequence,  this  mapping  is  performed  prior  to  the  conversion  to  pressure  so  that  the  SI 
correction  from  Eq.  (6)  can  be  performed.  Alignment  of  the  two  sets  of  raw  images  is  not  needed  in  this 
case. 
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Figure  3.  Basic  Processing  Sequences 

Dashed  lines  are  used  in  Fig.  3  to  indicate  which  macros  have  to  be  performed  for  each  of  the  eight 
cameras  individually.  This  grouping  of  macros  has  implications  for  parallel  processing,  which  is 
discussed  in  Section  2.6.3. 

2.5  BENCHMARKS 

Two  data  sets  have  been  selected  to  serve  as  benchmarks  for  timing  studies  (see  Table  1).  These  are 
small  subsets  from  actual  PSP  tests  on  scale  models  of  an  F-16C  fighter  jet  and  NASA’s  X-38  Crew 
Return  Vehicle  (see  Fig.  4).  Both  tests  were  conducted  in  the  16-ft  Transonic  Wind  Tunnel  (16T)  at 
AEDC,  with  eight  cameras  mounted  in  the  test  section,  as  shown  in  Fig.  5.  (Only  six  cameras  were  used 
for  the  actual  data  reduction  on  the  X-38.)  In  both  tests,  FIB-based  paints  were  used  (see  Table  1)  as  the 
basis  for  an  intensity-based  measurement,  while  excitation  of  the  PtTFPP  fluorophore  was  accomplished 
with  xenon  arc  lamps,  which  were  filtered  between  350  and  550  nm. 


Figure  4.  Examples  of  Processed  Data  for  the  Two  Benchmark  Cases 
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Figure  5.  Schematic  of  Wind  Tunnel  Setup  with  Positions  of  the  Eight  Cameras  Shown 

In  both  benchmarks,  four  images  are  used  per  camera  per  data  point  (wind  off  and  wind  on,  each  with 
associated  black  images).  Image  sizes  are  1024  by  1024  pixels,  with  16-bit  resolution  per  pixel.  Both 
benchmarks  contain  four  model  attitudes,  taken  from  an  alpha  sweep.  The  3D  geometry  files  for  both  tests 
are  multizone,  structured.  Plot  3D  grids  with  over  300,000  grid  points  each.  In  both  cases,  self¬ 
illumination  corrections  were  performed  at  about  1/10  of  this  resolution,  with  about  30,000  cells  per 
calculation.  Further  details  are  given  in  Table  1.  Figure  6  shows  an  example  of  the  effect  of  the  SI 
correction  on  the  X-38. 


Table  1.  Benchmark  Data  Sets 


F-16C  (tf960) 

X-38  (tf976) 

Number  of  Data  Points 

4(a  =  0,  8,  18,  26  deg) 

4  (a  =  4,  8,  16,  20  deg) 

Mach  Number 

0.8 

0.6 

Number  of  Cameras 

8 

6 

Paint 

PtTFPP  in  FIB7 

PtTFPP  in  UniFIB 

Images  per  Camera  per  Point 

4 

4 

Image  Size 

1024x1024x16  bit 

1024x1024x16  bit 

Number  of  Registration  Targets 

102 

36 

Number  of  Pressure  Taps 

400 

118 

Number  of  Zones  per  Grid 

379 

264 

Number  of  Points  on  3D  Grid 

344,269 

311,726 

Number  of  Elements  in  A-Matrix 

29,172 

31,344 

Size  of  A-Matrix  File 

230  MBytes 

129  MBytes 

2.6  SOFTWARE  AND  HARDWARE 
2.6.1  Software 

Since  1998,  PSP  data  processing  at  AEDC  has  been  accomplished  with  “Green  Boot,”  a  program  that 
was  developed  in  the  mid  1990s  by  McDonnell  Douglas  Aerospace  (MDA)  and  NASA  Ames.  Since  the 
program  was  acquired  by  AEDC,  NASA  Ames  and  AEDC  have  continued  development  of  a  government 
version  of  the  code,  while  Boeing  (upon  acquiring  MDA)  has  continued  development  of  a  proprietary 
version.  The  government  version  of  the  code  (currently.  Version  Gb.2.17a)  consists  of  about  92,000  lines 
of  code,  programmed  88  percent  in  C  and  12  percent  in  Fortran.  The  program  has  an 
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Figure  6.  Illustration  of  Effect  of  Self-Illumination  Correction  on  X-38  Tail  Section 

extensive  graphieal  user  interfaee  (GUI)  but  also  supports  seript-based  proeessing.  The  program  runs  on 
either  SGI  Irix  or  Linux  platforms.  The  program  strueture  is  essentially  sequential,  though  some 
modifieations  have  been  made  to  support  parallel  proeessing,  as  deseribed  in  Seetion  2.6.3. 

Table  2  lists  the  prineipal  support  files  for  the  Green  Boot  program.  The  SQL  database  file  keeps 
traek  of  all  raw  and  proeessed  images  assoeiated  with  a  partieular  test,  as  well  as  assoeiated  parameters 
sueh  as  tunnel  eonditions  and  mapping  eoeffieients.  The  targets  file  eontains  model  eoordinates,  normal 
veetors,  and  (when  applieable)  diameters  of  eaeh  of  the  targets  and  pressure  taps,  as  well  as  other  points 
(“fidueials”)  that  need  to  be  loeated  in  the  image  but  that  are  not  used  for  image  registration.  The  eamera 
file  stores  parameters  that  are  required  to  aeeomplish  automatie  image  registration.  The  setup  fde  eontains 
various  program  settings  and  maintains  a  list  of  most-reeent  eommands.  The  wind  tunnel  data  (WTD) 
files  report  tunnel  eonditions  and  aerodynamie  eoeffieients,  as  well  as  pressure  readings  from  pressure 
taps,  if  any,  to  be  used  for  in-situ  ealibrations.  The  3D  grid  file  defines  the  geometry  of  the  test  artiele. 
The  eomponent  file  defines  groupings  of  zones  for  load  ealeulations  on  parts  of  the  model  (e.g.,  left  upper 
wing,  left  lower  wing,  fuselage,  ete.).  The  A-matrix  file  eontains  the  preealeulated  SI  eoeffieients  from 
Eq.  (6).  Maero  files  define  proeessing  sequenees  sueh  as  those  in  Fig.  3.  The  log  file  and  diagnosties  files 
maintain  a  reeord  of  submitted  inputs  and  outputs  generated  by  the  program. 


Table  2.  Support  Files  for  Data  Processing 


Input  Output  Calibration 

SQL  Database  File 

Targets  File 

Camera  File 

Setup  File 

WTD  Files 

3D  Grid  File 

Component  File 

A-Matrix  File 

Maero  Files 

Fog  File 

Diagnosties  Files 
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Table  2  indicates  whether  the  files  are  used  for  input,  output,  or  both,  and  whether  the  files  are  used 
as  part  of  the  manual  pretest  calibration  process,  which  has  to  be  completed  before  fully  automated 
processing  can  be  accomplished  in  run  mode. 

2.6.2  Hardware 

The  Green  Boot  code  has  been  run  on  a  variety  of  Silicon  Graphics  workstations,  from  the  now- 
obsolete  Indigo  2  and  02,  through  the  Origin  2000,  the  Octane,  and  Octane  2,  with  processors  ranging 
from  an  1P22,  R4400  for  the  Indigo  2  to  a  dual  IP  30,  R14000  for  the  Octane  2,  and  at  processing  speeds 
ranging  from  195  MHz  for  the  Indigo  2  to  600  MHz  for  the  Octane  2. 

Until  recently,  all  PSP  data  acquisition  (along  with  a  large  part  of  the  data  processing)  was  performed 
on  an  SGI  Origin  2000,  equipped  with  1.3  GB  of  memory  and  with  eight  1P27,  RIOOOO  processors 
running  at  the  (now  slow)  speed  of  1 95  MHz. 

Most  recently,  the  Green  Boot  code  has  been  transitioned  to  a  Linux  cluster.  This  system  is 
configured  to  allow  simultaneous  data  acquisition  and  data  processing.  It  consists  of  a  set  of  front-end, 
dual  P4,  Xeon  processors  running  at  2.4  GHz  and  eight  dual-processor  slaves  (also  P4  Xeons),  each 
running  at  2.2  GHz.  Communication  between  the  master  and  the  slaves  is  through  an  HP  switch  with  a  1- 
Gbit  Ethernet  connection  to  the  front  end  and  1 00-Mbit  connections  to  each  of  the  nodes.  The  front  end 
has  2  GB  of  memory,  while  the  nodes  have  1  GB  each.  Both  the  master  and  the  slaves  run  Linux  Red  Hat 

7.3  as  the  operating  system. 

2.6.3  Implementation 

Various  configurations  of  the  software  have  been  employed,  ranging  from  purely  sequential 
processing  (i.e.,  one  camera  at  a  time)  to  parallel  processing,  as  shown  in  Fig.  7.  The  implementation 
shown  in  Fig.  7  is  the  latest  and  fastest,  especially  if  a  separate  processor  is  available  for  each  of  the 
processes  shown  (i.e.,  one  master  process  and  up  to  eight  slave  processes). 


Figure  7.  Parallel-Processing  Implementation  for  the  Processing  Schemes  from  Figure  2.3 

The  implementation  from  Fig.  7  applies  to  both  processing  sequences  from  Fig.  3.  In  each  case, 
macros  in  the  dashed  boxes  in  Fig.  3  are  performed  in  parallel  for  each  camera  individually.  Upon 
completion,  partially  processed  data  from  all  cameras  are  combined  and  processed  to  final  data  on  the  3D 
grid. 

Because  the  Green  Boot  code  is  not  intrinsically  parallel  in  nature,  communication  between  the 
master  process  and  the  slave  processes  in  Fig.  7  is  accomplished  via  pipes  and  sockets  (see  Section  3): 
The  master  sends  out  a  start  signal  by  writing  a  run  and  sequence  number  into  a  signal  file  that  is 
monitored  continuously  by  the  slave  processes.  At  this  point,  the  master  process  itself  becomes  a  slave, 
waiting  for  a  completion  signal  (as  well  as  partially  processed  data)  from  each  of  the  original  slave 
processes,  and  also  via  the  hard  disk. 
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The  newer  software  program  (the  so-called  BlackBox  code)  accomplishes  the  same  tasks  by  use  of 
the  Parallel  Virtual  Machine  (PVM)  paradigm.  This  newer  program  achieves  greater  efficiency  for  on-line 
processing  by  allowing  more  efficient  data  sharing  between  processors  and  by  implementing  only  that 
functionality  of  the  GreenBoot  code  that  is  essential  for  on-line  processing  per  se  (i.e.,  without 
implementing  the  substantial  requirements  for  GUI  support,  pretest  calibration  functions,  and  posttest 
interrogation  of  the  processed  data). 


2,7  TIMING  RESULTS 

Benchmark  results  for  the  data  from  Section  2.5  are  presented  in  Tables  3  and  4.  In  these  tables,  the 
designations  “2D”  and  “3D”  refer  to  the  processing  schemes  from  Fig.  3.  In  each  case,  the  highest 
possible  optimization  settings  were  used  at  compile  time,  i.e.,  the  “-03  -rl2000”  flags  on  the  SGI 
machines  and  the  “-0”  flag  on  the  Linux  machines. 


Table  3.  Cumulative  Single-Processor  Times  (in  Seconds)  for  X-38  Data  on  SGI  Octane  2 


Step 

DEFINE 

0.4 

LOAD 

4.6 

REGISTER 

14.5 

ALIGN 

2.6 

CONVERT 

7.6 

MAP 

11.3 

SI-CORRECT 

SAVE 

0.5 

TOTAL 

41.4 

2D  Processing 


0  5  10  15 


3D  Processing 


0  5  10  15  20  25 

49.2 


Table  4.  Processing  Times  per  Data  Point  in  Seconds 


SGI  Origin  2000 

(8  Proc,  195  MHz) 

SGI  Octane  2 

(2  Proc,  600  MHz) 

Linux  Cluster 

(9  Dual  Proc,  2,2  GHz) 

F-I6C 

X-38 

F-I6C 

X-38 

F-16C 

X-38 

2D,  Seq 

150 

103 

60 

39 

32 

22 

3D,  Seq 

193 

121 

74 

48 

41 

28 

2D,  Par 

27 

22 

30 

22 

7 

6 

3D,  Par 

41 

31 

38 

28 

10 

8 

Table  3  shows  average  processing  times  per  processing  step  for  the  X-38  dataset  on  the  SGI  Indigo  2, 
when  only  a  single  processor  was  used.  It  is  seen  that  the  most  time-consuming  processing  steps  are 
image  registration  and  mapping  data  to  the  3D  grid.  Mapping  times  in  the  3D  case  are  twice  those  in  the 
2D  case,  which  is  consistent  with  the  fact  that  twice  as  many  images  have  to  be  mapped  (i.e.,  wind-off 
and  wind-on  intensities,  as  opposed  to  pressure  values). 

Table  4  lists  total  processing  times  per  data  point  for  both  the  F-16C  and  X-38  datasets.  Times  are 
grouped  by  machine,  by  processing  type  (2D  vs.  3D — see  Fig.  3),  and  by  sequential  (“seq”)  vs.  parallel 
(“par”)  processing.  In  the  sequential  case,  only  one  processor  is  used  per  machine.  This  results  in  longer 
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times  for  the  F-16C  data  than  for  the  X-38  data,  because  of  the  use  of  eight  cameras  versus  six.  In  the 
parallel  case,  either  eight,  two,  or  nine  processors  are  used,  as  indicated  by  machine  type.  Times  for  the  F- 
16C  and  X-38  are  closer  in  this  case.  Most  importantly,  parallel  processing  times  on  the  fastest  machine 
(the  Linux  cluster)  are  well  below  15  s  per  data  point.  This  time  (15  s)  was  the  target  for  achieving  on-line 
processing,  in  the  sense  that  at  this  speed,  one  data  point  can  be  processed  in  full  while  the  next  data  point 
is  being  acquired. 

2.8  CONCLUDING  REMARKS 

With  the  implementation  of  parallel  processing  of  PSP  image  data  on  a  Linux  cluster,  a  significant 
milestone  was  achieved-namely,  the  ability  to  reduce  processing  times  sufficiently  that  one  data  point  can 
be  processed  in  full  while  the  next  data  point  is  being  acquired.  Further  improvements  have  been  made, 
particularly  the  PVM-based  program  that  enables  more  efficient  sharing  of  data  among  multiple 
processors  than  was  previously  feasible  with  the  Green  Boot  code. 

The  processing  schemes  employed  can  easily  be  extended  to  more  general  PSP  measurements-for 
example,  those  that  would  employ  a  third  (and  possibly  a  fourth)  gate  for  simultaneous  measurement  of 
pressure  and  temperature,^®'*^  or  for  a  temperature-corrected  PSP  measurement. 

Further  reductions  in  data-processing  times  (with  processing  times  moving  below  one  second  per  data 
point)  are  probably  not  required  for  the  study  of  steady-state  flows.  However,  it  is  conceivable  that  further 
speed  increases  might  become  desirable  for  the  study  of  unsteady  flow  phenomena  such  as  buffet,  flutter, 
limit-cycle  oscillations,  and  other  dynamic  pressure  phenomena.  Such  further  improvements  in  data- 
processing  times  would  probably  call  for  both  hardware  and  software  improvements  beyond  those 
described  here,  or  would  use  two-tier  processing,  with  approximate  data  being  generated  in  near-real 
time.  More  accurate  processing  would  then  be  performed  after  the  completion  of  the  wind  tunnel  test. 
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3.0  MESSAGE  PASSING  FOR  PARALLEL  PROCESSING  OF 
PRESSURE-SENSITIVE  PAINT  IMAGES^ 


3.1  ABSTRACT 

A  message-passing  scheme  is  described  that  allows  parallel  processing  of  pressure-sensitive  paint 
images  on  a  machine  with  multiple  processors  or  a  cluster  with  multiple  nodes.  The  scheme  implements 
the  use  of  forks  and  pipes  in  the  former  case  and  socket-based  TCP/IP  communications  in  the  latter.  The 
approach  demonstrates  how  multiple  copies  of  a  nonparallel  legacy  code  (in  this  case,  NASA’s  Green 
Boot  software)  can  be  made  to  run  in  parallel  in  either  a  parent-child  or  a  client-server  configuration. 
Results  are  presented  for  benchmark  data  from  wind  tunnel  tests  of  an  F-16C  fighter  jet  model  and 
NASA’s  X-38  Crew  Return  Vehicle. 


3.2  INTRODUCTION 

Pressure-sensitive  paint  (PSP)  has  established  itself  as  an  important  test  and  evaluation  tool  for 
measuring  full-field  pressure  distributions  on  test  articles  in  aerodynamic  test  facilities,  particularly 
transonic  wind  tunnels. '  At  the  AEDC,  these  pressure  distributions  are  obtained  by  processing  image  data 
from  up  to  eight  digital  cameras  in  the  test  section  of  the  wind  tunnel  (see  Fig.  5).  Final  data  are  presented 
to  the  customer  on  a  three-dimensional  (3D)  grid  of  the  test  article.  Processing  steps  include  automatic 
target  detection,  image  registration,  image  alignment,  reflected-light  correction,  conversion  of  signal 
ratios  to  pressure,  and  mapping  of  image  data  to  the  3D  grid.^^’^"* 

Data  processing  is  accomplished  with  Green  Boot,  a  code  that  was  developed  originally  by  NASA 
Ames  Research  Center  and  McDonnell  Douglas  Aerospace  (MDA).  The  code  consists  of  some  90,000 
lines  of  C  (with  some  Fortran),  offers  extensive  graphical  user  interface  (GUI)  support,  and  has  provisions 
for  script-based  processing.  It  has  always  been  possible  to  run  multiple  copies  of  the  code  as  separate 
processes.  However,  management  of  the  required  databases  and  script  files  has  proven  to  be  cumbersome 
and  an  impediment  to  achieving  automated  processing  in  near-real  time  (i.e.,  within  15  s  of  taking  data). 


3.3  PARALLELIZATION  APPROACH 


Figure  8  illustrates  how  data-processing  requirements  naturally  suggest  a  master-slave 
implementation  of  parallel  processing  PSP  image  data:  a  master  process  initiates  a  request  to  process  a 
data  point,  whereupon  up  to  eight  slave  processes  (one  for  each  camera)  perform  partial  processing  of 
image  data  from  a  particular  camera-the  result  being  processed  data  that  are  mapped  to  the  3D  grid.  The 
master  process  then  collects  this  mapped  data  from  the  slaves,  merges  the  data  onto  a  single  3D  grid,  and 
completes  processing  of  the  data.  These  steps  are  described  in  detail  in  Ref  24.  The  present  section 
focuses  on  the  message-passing  scheme  that  was  developed  to  implement  this  concept. 


Master  Process 


Initiate  Request  to 
Proeess  Data  Point 

Merge  Data  on 
Single  3D  Grid 


Complete  Proeessing 
and  Save  Result 


Slave  Processes 
(e.g.,  one  per  eamera) 


p 

'l 

-_ 

— 

1  , 

■_ 

— 

=1-1 

Per-Camera 
Processing: 
Data  on  3D 
Grid 

Figure  8.  Master-Slave  Approach  to  PSP  Data  Processing 


'^W.  Ruyten  and  W.  E.  Sisson,  “Message  Passing  for  Parallel  Processing  of  Pressure-Sensitive  Paint  Images,”  DoD 
HPCMP  Users  Group  Conference,  Williamsburg,  Virginia,  June  7-11,  2004.  (Also  published  by  IEEE  in  bound 
Proceedings,  ISBN  0-7695-2259-9,  pp.  308-312). 


19 


AEDC-TR-06-6 


Figure  9  lists  a  typical  macro  that  is  executed  by  the  master  process:  It  defines  the  names  of  the  to-be- 
processed  images  and  the  final  3D  data  file  on  the  basis  of  run  and  sequence  numbers  of  the  data  points.  It 
then  instructs  each  of  the  slave  processes  in  Fig.  8  to  execute  the  script  slave_macro  for  a  particular 
camera,  merges  the  returned  data  onto  the  3D  grid,  and  saves  the  final  result  in  Plot-3D  (or  other)  format. 
The  INSTRUCT  and  MERGE_3D  commands  are  built  into  the  Green  Boot  code,  whereas  _def  ine  and 
_save  are  script  files  that  are  defined  by  the  user  in  terms  of  native  Green  Boot  commands.  The  variable 
$P3D  is  set  by  the  script  _def  ine. 

In  Figure  9,  the  instruction  “INSTRUCT  $CAM  ...”  causes  the  script  “slave_macro  $I  $2 
$  CAM”  to  be  executed  on  the  slave  process  assigned  to  camera  $  CAM.  (Here,  $  I  and  $  2  represent  the  run 
and  sequence  numbers  of  the  data  point.)  Figure  10  lists  a  typical  realization  of  slave  macro.  It  loads  the 
reference  and  run  images  from  the  specified  camera;  finds  the  registration  targets;  registers  the  images  by 
determining  the  mapping  transformation  from  2D  image  space  to  3D  model  space;  aligns  the  two  images 
(to  account  for  small  rotations  and/or  shifts  between  the  two);  ratios  the  images;  converts  the  ratio  to  a 
pressure  or  pressure  coefficient  (based  on  a  calibration  equation);  maps  the  resulting  pressure  to  the  3D 
grid;  and  sends  this  information  back  to  the  master  process.  The  variables  $REF,  $RUN,  and  $P3DCAM 
in  Fig.  10  are  each  set  by  the  script  “_def  ine.”  All  scripts  prefaced  with  an  underscore  are  defined  by 
the  user. 

The  processing  sequences  from  Figs.  9  and  10  can  be  used  both  with  intensity-based  PSP  (in  which  a 
reference  image  is  obtained  with  the  tunnel  at  atmosphere)  and  with  lifetime-based  PSP,  in  which  both  the 
reference  and  the  run  images  are  obtained  at  the  run  condition.  Variations  on  these  scripts  are  possible. 
For  example,  in  order  to  perform  a  correction  for  reflected  fluorescent  light,  the  reference  and  run  images 
both  have  to  be  fully  mapped  to  the  3D  grid  before  conversion  to  pressure  can  take  place.  In  this  case,  the 
master  macro  performs  MERGE_3D  on  two  sets  of  images,  performs  the  reflected-light  correction,  and 
only  then  performs  the  _c  on  vert  macro  before  saving  the  final  result. 


define  $1  $2 

// 

$1,  $2  define  data  point 

foreach  CAM  (01 

02  .  . 

.  08) 

// 

Loop  over  all  cameras 

INSTRUCT  $CAM 

slave 

macro  $1  $2  $CAM 

// 

Process  $CAM  on  Slave  $CAM 

end 

// 

End  of  loop 

MERGE  3D  $P3D 

// 

Combine  data  on  3D  grid 

save  $P3D 

// 

Save  final  result 

Figure  9.  Typical  Master  Macro  Script 


define  $1  $2  $3 
Toad  $REF  $3 
Toad  $RUN  $3 
register  $REF 
register  $RUN 
"align  $RUN  $REF 
convert  $RUN  $REF 
map  $RUN  $P3DCAM 
"send  $P3DCAM 


//  $1,  $2,  $3  define  data  point,  camera 
//  Load  reference  image  (e.g.,  wind-off) 

//  Load  run  image  (e.g.,  wind-on) 

//  Find  targets  and  register  ref  image 
//  Find  targets  and  register  run  image 
//  Align  the  run  image  to  the  ref  image 
//  Ratio  images  and  convert  to  pressure 
//  Map  pressure  data  to  3D  grid  for  this  camera 
//  Send  result  back  to  master  process 


Figure  10.  Typical  Realization  of  slave_macro  in  Figure  3.2 


3.4  IMPLEMENTATION 

The  master-slave  communications  scheme  from  Fig.  8  has  been  implemented  in  two  ways.  In  the  first 
approach,  the  master  process  is  forked  repeatedly  before  X-window  support  is  requested,  and  write  and 
read  pipes  are  established  from  the  parent  process  (the  master)  to  the  resulting  child  processes  (the 
slaves).  In  the  second  approach,  slave  processes  are  started  up  in  server  mode  on  one  or  more  host 
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machines  at  preselected  ports,  and  the  master  process  acts  as  a  client  that  requests  a  socket-based  TCP/IP 
connection  to  each  slave.^^  In  both  cases,  a  set  of  fde  descriptors  is  established  on  each  process  (f  dw  [  ] 
for  writing,  f  dr  [  ]  for  reading),  as  illustrated  in  Fig.  11.  This  allows  the  program  to  communicate 
between  processes  by  performing  system-defined  WRITE  and  READ  calls  on  these  file  descriptors,  as 
illustrated  in  Fig.  12  for  a  set  of  data  contained  in  a  buffer,  buf .  The  while  loop  in  Fig.  12  ensures  that 
data  packets  are  received  in  full,  even  if  more  than  one  set  of  write/read  operations  is  required.  Typically, 
0.9  to  1.4  MB  of  data  are  returned  to  the  master  process  per  slave  process  per  data  point.  These  data  are 
buffered  in  packets  of  up  to  10,000  bytes,  a  size  that  appears  to  work  well  for  both  the  pipe  and  the  socket 
implementations. 


Master 

Slave  1 

fdr [0] 

fdw[0] 
fdr [0] 

f dw  [  0 ] 

fdw[l] 
fdr [1] 

Slave  2 

fdr [0] 
f dw  [  0 ] 

fdw[2] 

Slave  3 

fdr [0] 

fdr [2] 

etc. 

f dw [ 0 ] 

etc . 


Figure  11.  Implementation  of  Master-Slave  Communication  Through  File  Descriptors 


int  len  =  sizeof  (buf) ; 

// 

Length  of  data 

in  buf 

write (fdw [0] , &len, sizeof (len) ) ; 

// 

Write 

length  of  data 

write (fdw [0 ] , buf , len) ; 

// 

Write 

data 

int  len,  nr=0; 

// 

Initialize  num 

bytes  read 

read ( fdr [ s] , &len, sizeof (len) ) ; 

// 

Read 

length  of 

data  from  “s” 

while (nr ! =len ) 

// 

Loop 

until  all 

data  read 

nr+=read(fdr[s] ,&buf[nr] ,len-nr) ; 

// 

Read 

data  from 

slave  “s” 

Figure  12,  Interprocess  Communication  via  System-Level  write  and  read  Calls 


Figures  13  and  14  show  partial  source  code  for  the  establishment  of  the  fde  descriptors  in  the  two 
scenarios.  The  fork/pipe  implementation  in  Fig.  13  is  straightforward  and  requires  separate  pipes  for 
writing  and  reading.  The  socket  implementation  from  Fig.  14  is  somewhat  more  involved,  but  it  allows 
both  reading  and  writing  on  a  single  socket.  Byte  swapping,  if  required,  is  performed  by  the  master 
process  only.  (For  example,  in  Fig.  12,  the  variable  len  would  be  byte-swapped  prior  to  the  first  write 
operation.) 


for  (int  s=0;  s<N;  s++)  { 

// 

Request  N  slave  processes 

int  MS [2] , SM[2] 

,  pid; 

if  (pipe(MS)<0 

II  pipe(SM)<0)  error... 

// 

Create  the  pipes 

if  ( (pid=f ork ( ) ) <0 )  error... 

// 

Fork  the  main  process 

if  (pid==0)  { 

// 

Parent  process  (master) 

fdw [ s ] =MS [ 1 ] ; 

close (MS [0] )  ; 

// 

Set  file  descriptors 

fdr [s ] =SM [0 ] ; 

close (SM [ 1 ] ) ; 

continue; 

I 

else  if  (pid>0) 

{ 

// 

Child  process  (slave) 

fdw[0]=SM[l] ; 

close (SM[0] ) ; 

// 

Set  file  descriptors 

fdr[0]=MS[0] ; 

close (MS [ 1 ] ) ; 

(enter  slave 

} 

} 

mode . .  . ) 

// 

Indefinite  loop 

Figure  13,  Fork-Pipe  Implementation  with  Master  and  Slave  Processes 
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Start  slave  process  in  server  mode: 

struct  sockaddr  in  sock;  int  fd; 

(populate  sock:  AF  INET,  port  number,  INADDR  ANY) 

// 

Specify  port 

int  fd  tmp  =  socket (AF  INET, SOCK  STREAM, IPPROTO  TCP) ; 

// 

TCP-IP  socket 

if  (bind(fd  tmp, &sock, sizeof (sock) ) <0)  error... 

if  (listen (fd  tmp, MAXQUEUE ) <0 )  error... 

// 

Set  queue  size 

if  ( ( f d=accept ( f d  tmp, NULL, NULL) ) <0)  error... 

// 

Connect  client 

fdr [0] =fdw [0] =fd;  close (fd  tmp) ; 

// 

File  descriptors 

(enter  monitoring  mode...) 

Start  master  process  as  client  to  slave  s: 

struct  sockaddr  in  sock;  struct  hostent  *host; 
if  ( (host=gethostbyname (hostname) ) ==NULL)  error... 

// 

Specify  host 

(populate  sock:  AF  INET,  port  number,  host->h  addr) 

// 

Specify  port 

int  fd  =  socket (AF^INET, SOCK  STREAM, IPPROTO  TCP) ; 

// 

TCP-IP  socket 

if  (connect (fd, &sock, sizeof (sock) ) <0)  error... 

// 

Connect  server 

fdr [s] =fdw [s] =fd; 

(continue  with  main  program. . . ) 

// 

File  descriptors 

Figure  14,  TCP/IP  Implementation  with  Master  and  Slave  Processes 


The  fork/pipe  approach  is  particularly  suited  to  a  single  machine  with  multiple  processors,  such  as  an 
SGI  architecture.  It  has  the  advantage  that  the  resulting  slave  processes  are  tied  directly  to  the  master 
process,  so  that  separate  management  of  the  slave  processes  is  not  required.  The  socket  approach  is 
intended  for  use  on  a  cluster  of  machines  in  which  multiple  processors  communicate  via  an  Ethernet 
connection.  In  this  case,  one  must  exercise  care  to  ensure  that  the  slave  processes  are  running  and 
communicating  when  the  master  requests  that  a  data  point  be  processed.  In  either  case,  the  resulting 
scheme  is  much  more  robust  than  the  one  that  was  used  in  the  prototype  version  described  in  Ref  3, 
which  relied  on  the  use  of  instruction  fdes  that  were  written  to  and  read  from  a  common  file  system.  In 
particular,  conflicts  that  arose  when  one  process  was  trying  to  read  from  the  file  system  while  another  was 
in  the  process  of  writing  to  it  (or  vice  versa)  are  avoided  with  both  the  fork/pipe  and  socket  approaches, 
which  perform  buffering  and  synchronization  of  read  and  write  operations  implicitly.  Master  and  slave 
processes  still  access  a  common  database  for  tracking  information  related  to  the  raw  and  processed  image 
data.  This  is  a  commercially  available  SQL  database  designed  specifically  to  synchronize  access  by 
multiple  processes. 

To  facilitate  management  of  the  slave  processes  on  a  cluster,  a  Green-Boot-specific  daemon  process 
is  registered  on  each  node  of  the  cluster  under  the  super-daemon  inetd.^^  This  allows  a  Green  Boot  slave 
process  to  be  started,  polled,  and  shut  down  remotely  from  the  node  on  which  the  master  process  is 
executed  (typically,  the  front  end  of  the  cluster). 

Some  further  details  on  the  operational  modes  of  the  improved  Green  Boot  code  are  presented  in 
Section  3.6. 

3,5  RESULTS 

Correct  functioning  of  the  code  was  verified  on  several  platforms,  including  an  SGI  Octane  2  (with 
dual  IP30,  R14000  processors  running  at  600  MHz),  an  older  SGI  Origin  2000  (with  eight  IP27,  RIOOOO 
processors  running  at  195  MHz),  and  a  Linux  cluster  consisting  of  a  front  end  and  eight  nodes.  The  Linux 
cluster  (running  Red  Hat  Linux  7.3)  is  powered  by  dual  P4  Xeon  processors  on  both  the  front  end  and  the 
nodes,  running  at  2.4  GHz  and  2.2  GHz,  respectively.  Communication  between  the  front  end  and  the 
nodes  is  through  an  HP  switch  with  a  1 -Gbit/s  Ethernet  connection  to  the  front  end  and  100-Mbit/s 
connections  to  each  of  the  nodes.  Measured  transfer  rates  on  the  cluster  (with  one  master  process  and 
eight  slave  processes  on  the  nodes)  exceed  80  Mbit/s. 
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When  the  code  was  run  on  a  single  multiprocessor  machine,  no  significant  difference  in  total 
processing  time  was  found  between  the  fork/pipe  and  socket  implementations.  This  is  consistent  with  the 
fact  that  interprocess  communication  requires  far  less  CPU  time  than  does  the  actual  processing  of  the 
data.  As  expected,  best  performance  was  obtained  on  the  Linux  cluster,  which  has  the  highest  processing 
speed  among  the  platforms  tested.  Total  processing  times  per  data  point  for  two  benchmark  datasets 
(involving  an  F-16C  fighter  jet  and  NASA’s  X-38  Crew  Return  Vehicle)  are  detailed  in  Section  2.7  and  in 
Ref.  24  for  the  three  machines  quoted  above,  and  are  shown  in  Table  5,  for  both  “2D”  processing  (without 
a  reflected-light  correction)  and  “3D”  processing  (with  a  reflected-light  correction).  In  the  F-16C  case, 
image  data  from  eight  cameras  were  processed  on  a  master  process  supported  by  eight  slave  processes; 
for  the  X-38  case,  image  data  from  six  cameras  were  processed  on  a  master  process  supported  by  six 
slaves.  In  both  cases,  four  1024x1 024x1 6-bit  images  are  used  for  each  camera  (wind-off  and  wind-on 
images,  each  with  an  associated  black  image  that  is  subtracted  as  part  of  the  load  macro  in  Fig.  10),  and 
final  data  are  mapped  onto  a  3D  grid  with  more  than  300,000  grid  points.  Figure  4  shows  sample  results 
for  the  two  benchmarks,  with  color  representing  the  value  of  the  measured  pressure.  (Red  represents  high 
pressure;  blue  represents  low  pressure.) 

From  a  performance  perspective,  the  results  in  Table  5  are  particularly  significant  in  that  a  single  data 
point  can  be  processed  in  well  under  1 5  s  on  the  Linux  cluster.  This  was  the  AEDC  goal  for  achieving 
near-real-time  processing  of  PSP  test  data. 


Table  5.  Parallel  Processing  Times  per  Data  Point  in  Seconds 


SGI  Origin  2000 

(8  Proc,  195  MHz) 

SGI  Octane  2 

(2  Proc,  600  MHz) 

Linnx  Clnster 

(1+8  Proc,  2.2  GHz) 

F-16C  X-38 

F-16C  X-38 

F-16C  X-38 

2D  (w/o  refl) 

3D  (w/  refl) 

27  22 

41  31 

30  22 

38  28 

7  6 

10  8 

3.6  CONCLUSION 

It  has  been  successfully  demonstrated  that  it  is  possible  to  parallellize  an  intrinsically  nonparallel 
legacy  code,  particularly  one  that  lends  itself  to  a  master-slave  configuration  in  which  the  master  process 
initiates  a  request  to  process  a  data  point,  slave  processes  perform  partial  processing  of  the  data,  and  the 
master  process  completes  processing  of  the  data  after  combining  the  results  from  the  slaves. 
Communication  between  the  master  and  the  slaves  can  be  accomplished  either  across  pipes  by  repeatedly 
forking  the  main  process  (on  a  multiprocessor  machine)  or  (on  a  cluster)  by  establishing  socket-based 
TCP/IP  communications  between  master  and  slaves,  with  each  slave  running  as  a  server  process  on  a 
node,  and  the  master  process  running  as  a  single  client  on  the  front  end  of  the  cluster.  This  approach  to 
parallelization  may  be  applicable  to  data  processing  schemes  other  than  those  used  here  for  pressure- 
sensitive  paint. 

3.7  APPENDIX  3A:  GREEN  BOOT  RUN  MODES 

The  improved  Green  Boot  code  can  be  run  in  four  modes,  depending  on  the  presence  of  an  optional 
master-slave  argument  on  the  command  line: 

(1)  “$PATH/gb  config”:  Start  a  single  process  on  the  specified  configuration  files  (i.e.,  a  collection  of 
data  files  with  different  extensions,  each  of  which  has  “config”  as  the  root  name  of  the  file.  (This  is  the 
original  Green  Boot  mode.) 

(2)  “$PATH/gb  config  -forkN”:  Start  a  master  process  that  will  fork  N  slave  processes.  The  slaves  are 
terminated  upon  exiting  the  master  process. 
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(3)  “$PATH/gb  config  -portP  Start  (in  background  mode)  a  slave  process  that  will,  upon  being 
contacted  by  a  master  client  process,  establish  a  TCP/IP  connection  to  the  master  process  on  (slave)  port 
number  P.  This  command  can  be  performed  either  manually  on  the  machine  upon  which  the  slave  process 
is  to  be  run,  or  automatically,  as  the  slave  process  can  be  spawned  by  the  inetd  super  server,  when  the 
super  server  receives  a  remote  request  from  the  master  process.  The  slave  process  is  terminated  when  the 
server  receives  the  built-in  Green  Boot  command  “QUIT”  from  the  master  process. 

(4)  “$PATH/gb  config  -serv”:  Start  a  master  process  that  will  request  that  a  series  of  slave  processes 
be  spawned  remotely,  depending  on  the  contents  of  a  setup  file.  The  setup  file  contains,  for  each  camera, 
the  hostname  and  port  number  of  the  node  that  is  to  run  a  slave  process  for  that  camera  in  server  mode,  as 
well  as  the  paths  to  the  database  files  and  the  Green  Boot  executable  on  the  node.  When  the  master 
process  is  able  to  verify  that  all  the  slaves  are  in  listening  mode  (see  Fig.  14),  the  socket  connections  are 
established.  When  the  user  exits  the  master  process,  the  master  sends  a  QUIT  command  to  each  of  the 
slave  processes. 
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4  .0  AUTOMATIC  IMAGE  REGISTRATION  FOR  OPTICAL  TECHNIQUES  IN 
AERODYNAMIC  TEST  FACILITIES* 

4.1  OVERVIEW 

A  mathematical  model  is  described  for  automatic  registration  of  digital  images  of  sting-mounted  test 
articles  in  aerodynamic  test  facilities.  The  model  is  based  on  three  premises:  the  use  of  registration  targets 
on  the  test  article;  approximate  knowledge  of  the  attitude  angles  of  the  test  article  with  respect  to  the  test 
facility  at  each  measurement  point;  and  the  ability  to  calibrate  the  cameras  with  respect  to  the  test  facility. 
Registration  of  individual  images  proceeds  in  three  steps:  first,  initial  predictions  are  calculated  for  the 
image  coordinates  of  all  visible  registration  targets;  then  template -based  searches  are  performed  to  locate 
the  targets  to  subpixel  accuracy;  finally,  an  exact  mapping  between  2D  image  coordinates  and  3D  model 
coordinates  is  established.  The  math  model  includes  deflection  of  the  sting  under  air  loads.  Sample  results 
are  presented  for  several  windtunnel  tests. 

4.2  INTRODUCTION 

A  variety  of  optical  techniques  are  used  to  study  aerodynamic  flows  in  test  facilities  and  their  effects 
on  vehicle  performance.  These  techniques  include  pressure-  and  temperature-sensitive  paints  (PSP,  TSP), 
optical  model  deformation  measurement,  particle -imaging  velocimetry  (PIV),  Doppler  global  velocimetry 
(DGV),  laser-induced  fluorescence  (LIT),  and  others  (see  any  of  the  proceedings  of  the  International 
Congress  on  Instrumentation  in  Aerospace  Simulation  Facilities,  e.g..  Ref  13).  A  common  challenge  for 
these  optical  techniques  is  to  establish  the  correspondence  between  2D  image  coordinates  and  the 
corresponding  coordinates  in  3D  space.  For  example,  in  the  case  of  PSP,  it  is  customary  to  map  measured 
pressures  to  a  3D  grid  of  the  test  article.  * 

Photogrammetry  is  a  well-established  technique  used  to  determine  the  required  mapping 
transformations  between  2D  image  coordinates  and  3D  spatial  coordinates  and  is  based  on  the  use  of 
registration  targets. The  3D  coordinates  of  these  targets  are  measured  before  the  test,  for 
example,  using  a  coordinate-measuring  machine.  The  challenge,  then,  is  to  locate  these  targets  in  the 
hundreds  or  thousands  of  images  that  are  acquired  as  part  of  the  diagnostic  employed,  while  the  test 
article  is  moved  about  the  test  facility.  The  purpose  of  Section  4.0  is  to  describe  a  scheme  that  has  been 
developed  at  the  Arnold  Engineering  Development  Center  (AEDC)  to  automate  this  image  registration 
process  for  digital  images  of  sting-mounted  test  articles.  Parts  of  this  image  registration  scheme  have  been 
described  elsewhere,  including  the  use  of  photogrammetry,'^’^*  the  initial  predictions  of  target  location  in 
images,*^’'®  the  use  of  template -based  searches  to  determine  the  actual  target  locations  in  an  image  to 
subpixel  accuracy  (see  Section  5.0),'^’^^  and  the  resulting  capability  to  perform  on-line  processing  of 
luminescent-paint  data  (see  Section  2.0).^"'  The  purpose  of  Section  4.0  is  to  address  those  aspects  that 
pertain  to  image  registration  per  se  in  a  comprehensive  fashion,  and  to  document  several  aspects  of  this 
process  that  are  not  covered  in  prior  publications  or  in  other  sections  of  the  present  report.  These  are  the 
use  of  a  sting  deflection  model  to  deal  with  the  displacement  of  a  test  article  caused  by  bending  of  the 
sting  under  the  effect  of  air  loads  in  a  wind  tunnel  (Section  4.4);  the  implementation  of  the  template -based 
search  scheme  (Section  4.6);  aspects  of  the  resulting  image  registration  algorithm  (Section  4.8);  and 
selected  results  (Section  4.9).  First,  the  use  of  photogrammetry  is  reviewed  in  Section  4.3. 

4.3  REVIEW  OF  PHOTOGRAMMETRY 

In  its  most  basic  form,  imaging  by  a  camera  is  described  as  a  projective  transform,  whereby  the 
spatial  coordinates  x,*  of  a  target  i  give  rise  to  a  pair  of  image  coordinates  (m„  v,),  given  by  (see  Fig.  15): 


*W.  Ruyten,  “Automatic  Image  Registration  for  Optical  Techniques  in  Aerodynamic  Test  Facilities,”  AIAA  Paper 
2004-2400,  24th  AIAA  Aerodynamic  Measurement  Technology  and  Ground  Testing  Conference,  Portland,  Oregon, 
June  28-July  1,  2004  (outstanding  paper  award). 
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Here  Uc  and  Vc  denote  the  image  coordinates  of  the  optical  axis  of  the  camera,  fc  is  the  effective  focal 
length  of  the  camera,  and  (7„  Vi,  and  Wi  are  the  components  of  a  vector,  {7„that  is  given  by 

U,  ^R^Cx^-x^).  (8) 

The  vector  x^*  and  matrix  R^*  denote  the  position  and  orientation  of  the  camera  and  may  be 
expressed  as 

~  ’ '^c  ’  ®c )  •  (^) 

As  indicated,  the  orientation  matrix  R^*  depends  on  three  Euler  angles.  The  exact  choice  of  angles  is 
not  of  particular  importance  (see,  e.g..  Refs.  14  or  15)  so  long  as  the  rows  of  the  matrix  R^*  correspond  to 
an  orthonormal  set  of  vectors  {Uc*,  Vc*,  Wc*}  in  3D  space  that  defines  the  orientation  of  the  image  plane 
and  the  viewing  direction  of  the  camera.  As  indicated  in  Fig.  15,  the  convention  followed  here  is  that  the 
origin  of  the  image  coordinate  system  is  in  the  bottom  left  comer  of  the  image.  Alternatively,  the  image 
origin  might  be  defined  as  being  in  the  upper  left-hand  comer,  which  would  require  a  reversal  of  the 
directions  of  the  vectors  v^*  and  Wc*  (see  Fig.  1  of  Ref.  15). 


Tunnel 

Coord. 


Figure  15.  Schematic  of  Imaging  Process 

The  superscripts  “(0)”  in  Eq.  (7)  indicate  that  deviations  from  perfect  imaging  are  not  yet  included. 
For  example,  radial  lens  distortion  (quantified  by  a  factor  Kic)  produces  distorted  image  coordinates  that 
are  given  in  terms  of  the  undistorted  coordinates  from  Eq.  (7)  by 

Ui  =  mP  -  8ui  ,  Vj  =  vP  -  5^;  ,  (10) 

where 

(Ali  +  Ali ) ,  5^,-  -  Ki^A^i  (A^ +  A^ • ) ,  (1 1) 

and 

^ui  =  >  Ki  =  .  (12) 

Higher  order  distortion  terms  can  be  included  as  well.'"^ 
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The  reason  for  the  use  of  asterisks  in  Eqs.  (8)  and  (9)  is  to  distinguish  tunnel  coordinates  {with 
asterisks)  from  the  corresponding  model  coordinates,  denoted  without  asterisks.  Tunnel  coordinates  are 
measured  with  respect  to  the  test  facility;  model  coordinates  are  measured  with  respect  to  the  test  article. 
Model  coordinates  of  registration  targets,  denoted  x„  are  initially  known  only  with  respect  to  the  model 
coordinate  system.  In  Section  4.4,  the  task  is  undertaken  to  estimate  the  corresponding  tunnel  coordinates 
X,*  as  the  test  article  is  moved  about  the  test  facility.  Assuming  that  both  the  exterior  parameters  (xf,  yf, 
A*  it^e*  x^c*,  Wc*)  and  the  interior  parameters  (fc,  Uc,  Vc,  Kic)  of  the  camera  have  been  calibrated,  it  is  then 
possible  to  estimate  the  resulting  image  coordinates  (m„  v,)  of  visible  registration  targets  (Section  4.5);  to 
find  the  exact  locations  of  the  targets  using  template-based  correlation  matching  (Section  4.6);  and  to 
complete  the  image  registration  process  (Sections  4.7  and  4.8). 


4,4  STING  MODEL 

Neglecting  effects  of  model  deformation,  the  relation  between  the  model  coordinates  x,  of  a 
registration  target  and  its  corresponding  tunnel  coordinates  x,*  may  be  written  as 

^  ^  ^ 

Xf  =  X^  +  ’  (13) 


where  x^*  and  denote  the  position  and  orientation  of  the  test  article  at  an  attitude  setting  k,  in  terms  of 
three  position  values,  x**,  yk*,  Zk*  and  three  Euler  angles,  a**,  P^*,  (j)**  according  to 


*  /  *  *  *  X 


-^ki^k’^k’h) ■ 


(14) 


As  is  the  case  in  Eq.  (9),  different  definitions  of  the  pitch,  yaw,  and  roll  angles  (i.e.,  a^*,  P**,  and  (])*:*) 
are  possible  (see,  e.g..  Ref  30). 

For  a  pitch-roll  mechanism  such  as  that  in  AEDC’s  transonic  wind  tunnel,  the  calculation  of  x^-*  and 
Rk*  is  described  in  Ref  15  for  the  case  of  no-sting  deflection  under  air  loads.  To  incorporate  such  sting 
deflection,  the  assumption  is  made  that  the  sting  deforms  into  a  circular  arc  in  the  pitch  plane  of  the 
model.  This  situation  is  depicted  (in  exaggerated  fashion)  in  Fig.  1 6  for  the  case  of  zero  model  yaw,  zero 
model  roll,  and  zero  prebend  angle  (i.e.,  the  x-axis  of  the  model  is  aligned  with  the  sting).  In  this  case,  the 
sting  deflection  angle  is  given  by  Aas  =  a  -  as,  where  a  and  as  are  the  pitch  angles  of  the  model  and  the 
sting,  respectively,  as  reported  by  the  wind-tunnel  data  system.  Based  on  the  assumption  of  a  circular  arc, 
the  radius  of  curvature  of  the  sting  is  given  by  R  =  L/Aas,  where  L  is  the  effective  length  of  the  sting  over 
which  bending  takes  place.  The  displacements  Axs  and  Azs  of  the  model  resulting  from  sting  deflection, 
measured  along  and  perpendicular  to  the  sting,  respectively  (see  Fig.  16),  are  given  by 

Axg  =T-RsinA«s  =-^frA«s  ,  (15) 

and 

Azg  =R -RcosA«s  =yTA«§  .  (16) 


The  small-angle  approximations  assume  that  Aas  is  expressed  in  radians,  and  these  approximations 
are  justified  because  Aas  is  typically  less  than  3  deg  in  magnitude. 

Figure  17  illustrates  the  resulting  displacements,  Ax*^*  and  Az/*,  measured  along  the  tunnel  x*-  and 
z*-axes  for  a  PSP  test  of  an  F-16C  model  at  0.9  Mach  number.^  Solid  symbols  represent  model  position 
measurements  based  on  fully  processed  image  data  from  eight  cameras  (see  Ref  8). 

The  continuous  curves  in  Fig.  17  are  curve  fits  based  on  the  small-angle  approximations  from  Eqs. 
(15)  and  (16),  giving  an  effective  sting  length,  L,  of  63  ±  1  in.  This,  indeed,  is  the  approximate  length  of 
the  sting  in  this  particular  wind  tunnel  test. 
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C  =  Wind-On  Position  of  A 

D  =  Effective  Arc  Center 

EpEj  =  Origin  of  Model  Coord.  System 

Figure  16.  Schematic  of  Sting-Deflecting  Model 


Figure  17.  Sting  Deflection  of  F-16C  Model  Under  Air  Loads 


More  generally,  the  effects  of  mounting  the  model  at  a  prebend  angle  ao  and  the  effects  of  rolling  the 
sting  by  an  amount  (|)s  must  be  included.  So  far  it  appears  that,  as  far  as  the  initial  prediction  of  the 
resulting  image  coordinates  of  registration  targets  is  concerned,  the  effects  of  out-of-plane  bending  caused 
by  sideslip  forces  can  be  ignored,  at  least  up  to  yaw  angles  of  8  deg  relative  to  the  wind  vector.  The  sting 
bend  angle  Aas  can  then  be  estimated  according  to 

Aug  =  a  -  Uq  -  tan~^(cos^zis  tanug).  (17) 
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If  the  origin  of  the  tunnel  coordinate  system  is  chosen  to  be  at  the  wind-off  pitch-roll  center  of  the 
sting  mechanism  (point  A  in  Fig.  16),  the  resulting  position  vector  x*:*  and  model  orientation  matrix  R^* 
from  Eq.  (13)  are  given  by 

=R^(as)R;c(^s)[Axs  +R>,(ao  +  Aas)Xm]  (18) 

and 

R/t  =R>'(<^s)Rx(<^s)R>'(‘^o  +'^ots)Rm  •  (19) 

Here  R,,  and  R^  are  pitch  and  roll  matrices  about  y-  and  x-axes,  respectively,  by  the  angles  indicated. 
The  displacement  vector,  Axs,  is  given  in  terms  of  the  components  from  Eqs.  (15)  and  (16)  by 

Axs  =-(Axs,0,Azs).  (20) 

The  minus  sign  in  Eq.  (20)  is  required  because  of  the  particular  choice  of  tunnel  coordinates,  with  an 
x*-axis  pointing  upstream  and  a  z*-axis  pointing  downstream  (see  Fig.  16). 

Finally,  x^  and  R^  in  Eqs.  (18)  and  (19)  describe  the  position  and  orientation  of  the  test  article  at  the 
wind-off  reference  condition,  i.e.,  at  zero  pitch,  yaw,  and  roll.  Typically,  it  is  adequate  to  express  x„  and 
Rm  as 

Xm  =  (^m  A  z^),  R^  =  D(-l,-tl,-l) ,  (21) 

where  and  are  the  horizontal  and  vertical  positions  of  the  origin  of  the  model  coordinate  system  and 
D  is  a  diagonal  matrix  with  elements  as  indicated.  As  is  the  case  in  Eq.  (20),  the  minus  signs  in  Eq.  (21) 
are  a  result  of  the  choice  of  coordinate  system  in  Fig.  16. 

4,5  RESULTING  PREDICTION  MODEL 

In  combination,  Eqs.  (7)  through  (21)  constitute  a  complete  set  of  equations  for  the  prediction  of  a  set 
of  image  coordinates  (m„  v,)  from  a  set  of  model  coordinates  x„  based  on  four  sets  of  parameters:  a  set  of 
exterior  camera  parameters  p^*  =  jc*  A*  ‘t'c*  Re*  ®c*);  ^  set  of  interior  camera  parameters  =  {fc, 
Uc,  Vo  Kjc};  a  set  of  variable  sting  angles  a^-  =  {a,  as,  it)s};  and  a  set  of  fixed  tunnel  parameters  t  =  {ao,  x^, 
Zm,  L}.  Symbolically,  the  resulting  transformation  may  be  written  in  terms  of  a  mapping  function,  F,  as 

(Ui  (22) 

Of  course,  the  mere  fact  that  a  set  of  image  coordinates  (ui,  v,)  can  be  calculated  for  a  target  i  does  not 
guarantee  that  the  target  will,  in  fact,  be  visible  in  the  image.  Three  additional  conditions  must  be  met  to 
ensure  such  visibility: 

(1)  The  predicted  image  coordinates  must  fall  within  the  bounds  of  the  image. 

(2)  For  the  usual  case  of  an  opaque  model,  the  surface  normal  of  the  target  must  point  toward  the 
camera. 

(3)  The  target  must  not  be  occluded  by  other  parts  of  the  test  article  or  tunnel  mounting  structure. 

The  first  of  these  conditions  can  easily  be  verified  based  on  the  known  size  of  the  image,  i.e.,  by 
requiring  that 

\  +  s  KUj  <Nii  -s  ,  \  +  s  KVj  <Ny  -s .  (23) 

Here  (1,  1)  and  (A„,  Ay)  represent  the  coordinates  of  the  bottom  left  and  upper  right  comers  of  the 
image.  The  skirt  width,  s  (typically,  a  few  image  pixels),  is  applied  to  eliminate  from  consideration  targets 
that  may,  in  fact,  be  wholly  or  partially  outside  the  image. 
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The  second  visibility  condition  can  also  be  verified  easily  if  a  target  normal,  n„  is  specified  for  each 
target,  i.  If  n,  is  expressed  in  model  coordinates,  then  Ri:*n,  is  the  orientation  of  the  target  in  tunnel  space. 
The  second  visibility  condition  can  thus  be  expressed  as 


<0. 


(24) 


The  third  visibility  condition  (that  a  target  not  be  occluded)  is  more  difficult  to  ascertain  for  a 
complex  3D  test  article,  but  it  can  nevertheless  be  checked  if  a  detailed  3D  geometry  grid  of  the  test 
article  is  available.  (This  is  typically  the  case  in  luminescent  paint  measurements.)  In  particular, 
reinterpret  Eqs.  (7)  through  (21)  as  applying  to  grid  points  g,  as  opposed  to  target  points,  i.  It  is  then 
possible  to  compute  a  set  of  values  {Ug,  Vg,  -Wg)  for  each  grid  point  g,  where  Ug  and  Vg  are  projected 
image  coordinates  from  Eq.  (10)  and  -  Wg  is  the  negative  of  the  third  component  of  the  vector  U  from  Eq. 
(8).  This  value,  -Wg,  may  be  interpreted  as  the  distance  from  a  grid  point  Xg*  to  the  focal  plane  of  the 
camera. 

The  values  {Ug,  Vg,  -Wg)  can  be  arranged  in  a  depth  buffer  (usually  called  a  z-buffer^'),  in  which  each 
pixel  ip,  q)  holds  a  value,  Wbuf,  that  represents  the  distance  from  the  camera  to  the  part  of  the  test  article 
that  is  viewed  by  that  pixel.  If  such  a  depth  buffer,  Whuf  (p,  q),  can  be  computed,  the  third  visibility 
condition  for  a  target  i  can  be  written  as 

-  <  ^buf  +  ^  (25) 

where  (pi,  q,)  are  the  projected  coordinates  of  target  i  rounded  to  the  nearest  pixel  value,  -If,  is  the 
distance  value  for  target  i,  and  Sw  is  a  tolerance  parameter  that  prevents  the  rejection  of  legitimate  targets 
because  of  small  discrepancies  between  grid  coordinates  and  target  coordinates.  If  the  depth  buffer  is 
calculated  from  a  closed  surface  grid,  Eq.  (25)  implies  Eq.  (24). 


4.6  TEMPLATE-BASED  TARGET  SEARCH 

Using  the  model  from  Sections  4.3,  4.4,  and  4.5,  it  is  possible  to  predict  the  image  coordinates  (m„  v,) 
of  all  registration  targets  that  should  be  visible  in  an  image.  The  next  step  is  to  find  the  actual  target 
locations  within  some  search  radius  from  the  initial  predictions  and  to  determine  the  exact  locations  of  the 
targets  to  subpixel  accuracy.  Two  techniques  are  widely  used  for  such  searches:  blob-finding  and 
template-based  correlation  matching.  Blob  finders  typically  do  not  take  into  consideration  the  shape  of  the 
registration  target.  Instead,  they  isolate  a  contiguous  set  of  pixels  that  contrasts  strongly  with  the 
background.  The  center  of  the  blob  is  calculated  by  centroiding,  i.e.,  by  performing  an  intensity-weighted 
average  of  individual  pixel  locations  in  the  blob,^^  following  subtraction  of  the  background  signal. 

By  contrast,  template-based  techniques  employ  a  search  template  (see  Fig.  1 8)  that  matches  the  shape, 
size,  and  orientation  of  the  desired  image  feature.  This  allows  for  more  robust  searching  than  is  possible 
with  a  blob  finder  and  reduces  the  probability  that  legitimate  image  features  will  go  undetected  and  that 
unwanted  artifacts  will  be  mistakenly  identified  as  targets. 

A  typical  search  template  consists  of  an  array  of  numbers  between  0  and  1  (see  Fig.  1 9)  in  which  each 
value  represents  the  fraction  of  the  corresponding  image  pixel  that  falls  inside  the  projected  outline  of  the 
desired  registration  target.  This  template  is  stepped  across  the  image  in  one-pixel  increments  to  find  a 
point  of  maximum  correlation.^^  Subpixel  centering  can  be  achieved  by  subjecting  the  resulting  array  of 
correlation  values,  in  turn,  to  a  centering  determination  (for  example,  by  sampled  correlation 
interpolation).^"^ 

A  new  variation  of  template-based  correlation  matching  is  described  in  Refs.  17  and  29,  in  which  the 
conventional  search  template  (such  that  as  in  Fig.  19)  is  augmented  by  two  additional  arrays  (denoted 
here  by  /„  and  — see  Fig.  20),  which  represent  the  derivatives  of  the  occupied  pixel  fractions  fo  with 
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A  =  Actual  Center  +  Outline 
B  =  Predicted  Center  +  Outline 

Figure  18.  Illustration  of  Use  of  Search  Template 
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Figure  19.  Example  of  Conventional  Search  Template 
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Figure  20.  Augmented  Templates  Corresponding  to  Figure  4.5 


respect  to  subpixel  shifts  of  the  search  template  along  the  image  axes.  It  is  shown  in  Ref  17  that 
incorporation  of  the  derivative  templates  allows  for  more  accurate  subpixel  centering  of  registration 
targets  than  is  possible  with  either  blob  finding  or  conventional  template  matching  and  that  augmentation 
of  the  search  template  produces  higher  correlation  values  than  does  the  nonaugmented  case. 

Below,  the  implementation  of  both  conventional  and  augmented  template  matching  is  reviewed  for 
use  with  image  registration. 
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4,6.1  Calculation  of  the  Search  Template 


Let  {x,y}  denote  a  set  of  n  points  (with  j  =  1,  . . n)  on  the  perimeter  of  a  registration  target  i  that  has 
passed  the  visibility  tests  from  Eqs.  (23)  through  (25).  In  particular,  for  a  solid,  circular  target  with  radius 
r,  and  position  x,  in  model  space,  the  set  of  points  {xy}  may  be  calculated  as 


1;  COS 

(  271/A 

+  m^-  sin 

r  27i/-y 

1  n  ') 

[  n  )_ 

7  . 


(26) 


Here  1,  and  m,  are  unit  tangents  that  are  calculated  from  the  surface  normal  n„  with  n,  =  1,  x  m„  and  r, ' 
is  essentially  the  target’s  radius  r„  with  a  small  correction  that  renders  the  area  of  the  resulting  n-sided 
polygon  equal  to  that  of  the  circular  target: 


2n/  n 
sin(27i  /  n) 


(27) 


The  assumption  is  made  that  the  surface  is  flat  on  the  scale  of  the  target.  Using  Eq.  (22),  a 
corresponding  set  of  image  coordinates  {uy,  Vy)  may  be  calculated.  By  interpolation,  this  set  of  projected 
perimeter  points  can  be  expanded  into  a  A^-sided  polygon  (with  K  >  n),  such  that  each  line  segment  is 
contained  inside  a  single  image  pixel  (see  Fig.  21).  Eet  {uk-i,  Vk-\)  and  {uk  v^)  be  the  end  points  of  the  k-th 
line  segment  for  some  target  i,  with  (mq,  vq)  =  («&  va:),  and  let  {pk  qk)  denote  the  indices  of  the  pixel  in 
which  the  A:-th  line  segment  is  contained.  The  partial  contributions  of  line  segment  k  to  the  occupied  pixel 
fractions  fo(pk  q)  of  a  column  of  pixels  that  extends  downward  from  the  pixel  {pk  qk)  is  then  given  by 


A/oA;  iPk^<l) 


(Uk-\  -  Uk) 


Vk  +  Vk-\ 


■<lk 


Uk-l  Uk  1  ^min  —  q  ^  qk 


(28) 


•  :  Original  projected  points 
O :  Interpolated  points  on  pixel  boundaries 


Figure  21.  Calculation  of  Template  from  Polygonal  Curve 
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Here  is  the  lowest  value  of  qk  in  the  template.  Individual  contributions  can  be  either  positive 
or  negative.  However,  when  contributions  from  all  line  segments  k  are  added,  the  resulting  template 
values  will  be  in  the  range  0  <fo<  1,  provided  that  the  projected  outline  of  the  template  is  traversed  in  a 
counter-clockwise  direction. 

To  calculate  the  derivative  templates  that  are  required  for  augmented  template  matching,  observe  that 
the  change  in  the  value  Afok  for  a  pixel  (pt,  qk)  caused  by  an  infinitesimal  shift  5m  of  the  template  along  the 
M-axis  is  equal  to  (v^  -v^-i)  5m;  observe  also  that  the  change  in  Afok  caused  by  an  infinitesimal  shift  5v 
along  the  v-axis  is  equal  to  (uk-i-Uk)  5v.  The  partial  contributions  of  line  segment  k  to  the  derivative 
templates/,  and /,  can  thus  be  calculated  as 

¥uk iPk  ^<lk)^Vk-  Vk_i  ,  iPk ^qk)^  Uk_i  -  Uk  .  (29) 

When  the  partial  contributions  Afuk  and  Af^k  of  all  k  line  segments  are  summed,  the  resulting  arrays  of 
values  of /,  and /,  are  in  the  range  from  -1  to  +1  for  each  pixel  (see  Fig.  20).  Inclusion  of  a  skirt  with  zero 
values  of  /o,/,,  and  /,  along  the  edges  of  the  template  (as  shown  in  Figs.  19  and  20)  has  been  found  to 
improve  the  performance  of  the  template  search. 

4,6.2  Conventional  Template  Matching 

In  conventional  template  matching,  a  perfect  match  between  an  array  of  template  values  fo  and  an 
image  feature  (represented  by  an  array  of  intensity  values  g)  is  obtained  if  there  exists  a  set  of  integer 
shifts  (Ap,  Aq)  that  renders  the  image  values  g(p+Ap,  ^-l-A^)  exactly  equal  to  the  scaled  template  values  b 
+  cfoip,  q),  where  b  is  the  image  background  value  outside  the  template,  c  is  a  contrast  term  (which  can  be 
positive  or  negative),  and  p  and  q  are  pixel  indices  in  the  template.  More  generally,  a  measure  of  the 
match  between  the  image  feature  and  the  template  is  provided  by  the  normalized  correlation 

<kfg{Ap,Aq)=  ^  (with  0  <  Jjg  <1),  (30) 

where 


Sff^YXfQ(p,q)-{h)Y, 

(31) 

PS 

S fg  (Ap,  Aq)  =  X  l/o  q)  -  {fo )]  [sip  +  Ap,q  +  Aq)  -  (g)] , 

(32) 

PS 

Sgg  (Ap,  Aq)  =  Y.  [sip  +  ^P,q  +  ^q)  -  (g)P  • 

(33) 

PS 

{fo)  =  Y.foip^q)  /Zi’ 

(34) 

PS  /  PS 

{giAp,Aq))=Y.SiP  +  ^P,q  +  ^q)  /Zl- 

PS  /  PS 

(35) 

The  sums  in  Eqs.  (31)  through  (35)  extend  over  all  pixels  (p,  q)  in  the  template  and  hence  over  the 
equivalent  region  of  the  image,  as  determined  by  the  shifts  Ap  and  Aq. 

Perfect  correlation  between  the  template  and  the  image  is  signified  by  a  value  Jfg  =  1  for  some  set  of 
integer  shifts  (Ap,  Aq).  In  practice,  the  template  is  stepped  across  the  image  by  varying  the  integer  shifts 
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and  until  a  point  of  maximum  and  sufficient  correlation  is  found  (i.e.,  a  value  Jfg  above  some 
threshold  J^m,  with  Jmin<  1  )• 

4.6.3  Augmented  Template  Matching 

Assume  that  a  conventional  template  match,  as  described  in  Section  4.6.2,  has  been  performed  in  the 
vicinity  of  a  predicted  target  location  (m„  v,),  resulting  in  the  tentative  identification  of  the  desired 
registration  target  at  a  location  (M,+Ap,  v,+A^)  for  some  set  of  integer  shifts  Ap  and  A^.  Now  imagine  that 
the  template  is  shifted  by  a  further,  subpixel  set  of  shifts  (Am,  Av).  To  first  order,  this  is  equivalent  to 
replacing  the  occupied  pixel  fractions  f^ip,  q)  with  the  values  foip,  q)  +  fu(p,  q)Au  +  f^ip,  q)Av  ,  where  fu 
and /v  are  the  derivative  terms  computed  in  Section  4.6.1  and  Am  and  Av  are  the  to-be-determined  subpixel 
shifts.  Substitution  of  this  first-order  Taylor  expansion  for  the  template  into  Eqs.  (30)  through  (35)  leads 
to  an  optimization  problem  with  respect  to  the  unknown  subpixel  shifts  Am  and  Av,  namely,  to  maximize 
the  resulting  value  of  J/g.  Seven  additional  correlation  sums  arise  (see  Fig.  22)  and  are  denoted  by  Suu,  S^v, 
Syy,  Sfu,  Sfy,  Sug{Ap,  Aq)  and  5'„g(Ap,  Aq).  These  sums  are  defined  analogously  to  Eqs.  (30)  through  (35) 
with,  in  this  case,  {/'„)  =  (/()  =  0. 


Figure  22.  Correlation  Sums  Arising  in  Augmented  Template  Search 


For  the  special  case  in  which  the  derivative  templates  fu  and  f  are  orthogonal  to  the  occupied  pixel 
fractions  f  (i.e.,  the  case  in  which  Sp  =  Sp  =  0),  the  solution  of  the  resulting  optimization  problem  is 
obtained  for  subpixel  shifts  Am  and  Av  given  by  (see  Ref  17): 


Au" 

(  S’ 

^uu 

S’  ^ 

^uv 
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^ug 

Av, 
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'^vv  J 

[^^gj 

(36) 


In  terms  of  these  subpixel  shifts,  the  best  estimate  for  the  center  of  the  registration  target  is  now  given 
by  (M,-l-Ap-l-AM,  v,-l-Aq'-l-Av).  The  associated  augmented  normalized  correlation  value,  Jfg*,  is  given  in  terms 
of  the  conventional  value  Jfg,  from  Eq.  (30),  by 

fg  ^ fg  T  ^  J fg  ■  (32) 

^ff^gg 
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For  circular  registration  targets,  a  simple  way  to  ensure  that  the  condition  Sfu  =  5)^,  =  0  is  satisfied  is  to 
shift,  in  Section  4.6.1,  the  polygonal  outline  of  the  projected  target  such  that  its  center  is  located  at  either 
an  integral  or  a  half-integral  position  in  both  the  u  and  v  directions.  More  generally,  Eqs.  (36)  and  (37) 
may  be  replaced  with  somewhat  more  involved  formulas  (quadratic  in  Am  and  Av)  that  include  the 
nonzero  sums  Sfu  andAyi,  (see  Section  3.1  of  Ref  17). 

By  iterating  the  calculation  of  the  search  template  (this  time,  without  shifting  the  projected  target 
center  to  an  integral  or  a  half-integral  location)  and  recomputing  the  subpixel  shifts  Am  and  Av,  even  more 
accurate  centering  of  the  template  with  respect  to  the  image  is  possible  (see  Ref  17).  Also,  Ref  29  (see 
Section  5)  addresses  the  extension  of  the  augmented  template  technique  to  poorly  rendered  image 
features,  for  which  the  implicit  assumption  of  uniform  values  of  background  and  contrast  breaks  down.  In 
Refs.  17  and  29  it  is  also  shown  how  the  precision  of  the  augmented  centering  determination  can  be 
estimated  from  the  achieved  correlation  values. 

4.7  COMPLETION  OF  THE  IMAGE  REGISTRATION  TASK 

With  all  (or  nearly  all)  of  the  registration  targets  found  that  were  deemed  visible  in  the  image,  the  task 
remains  to  complete  the  image  registration  process.  The  exact  form  in  which  this  is  to  be  accomplished 
depends  on  the  application.  For  example: 

(1)  The  discrepancy  between  the  predicted  and  measured  image  coordinates  (m„  v,)  is  attributed  to  the 
position  vector  x^*  and  orientation  matrix  R^^*  from  Eq.  (14)  being  different  from  the  prediction  from  Eqs. 
(15)  through  (21).  This  calls  for  a  nonlinear  least-squares  fit  of  the  measured  image  coordinates  (m„  v,) 
using  Eqs.  (7)  through  (14),  with  Xk*,  yt*,  Zk*,  a**,  P**,  and  (])*:*  as  fit  parameters.  This  produces  an 
updated  mapping  transformation  between  2D  image  and  3D  model  coordinates,  as  well  as  a  measurement 
of  the  actual  location  and  attitude  of  the  test  article.  A  useful  measure  of  the  agreement  between  measured 
(“exp”)  and  fitted  (“fif ’)  image  coordinates  is  provided  by  the  root-mean-square  difference,  a,  between 
the  two,  defined  as 


a  = 

[-1 

2 

+ 

y(exp)  _  .^,(fit) 

ni 

hvl 

,  1  1  , 

,  1  1  , 

1 1 

where  N  is  the  number  of  registration  targets  in  the  image. 

(2)  To  improve  the  accuracy  of  an  optical  model  attitude  determination,  data  from  two  or  more 
calibrated  cameras  can  be  combined  in  a  single  least-squares  fit,  defined,  again,  in  terms  of  a  set  of  to-be- 
determined  parameters  Xk*,  yt*,  zt*,  a**,  P**,  and  This  approach,  including  the  estimation  of  the 
uncertainties  of  the  resulting  position  and  attitude  determination,  is  discussed  in  more  detail  in  Ref.  15. 

(3)  At  the  outset  of  a  wind  tunnel  test,  it  is  necessary  to  calibrate  the  exterior  and  interior  parameters 
of  the  cameras,  as  well  as  some  of  the  parameters  that  are  used  in  the  sting  model  from  Section  4.3,  most 
notably  the  offsets  Xm  and  in  Eqs.  (21)  and  the  effective  sting  length,  L,  in  Eqs.  (15)  and  (16).  See,  e.g.. 
Refs.  14  and  15.  During  this  calibration,  it  typically  is  necessary  for  the  user  to  select  a  number  of 
registration  targets  manually  when  setting  images. 

(4)  A  high  rms  error  in  Eq.  (38)  may  indicate  that  the  assumption  of  rigid-body  motion  is  invalid  and 
that  some  amount  of  model  deformation  occurs.  Nevertheless,  to  achieve  an  accurate  mapping  between 
2D  image  coordinates  and  3D  spatial  coordinates,  additional  parameters  must  be  introduced.  A  simple 
approach  (though  without  a  clear  theoretical  justification)  is  to  replace  the  projective  transform  equations 
from  Eqs.  (7)  and  (8)  with  a  direct  linear  transform  of  the  form^"'’'^’^° 

^(0)  _  L^Xj  +  L2Ti  +  L^Zj  +  L4 

'  LgX;  +  LjoTi  + +  1 


^(0)  ^  L^Xj  +  L(^yi  +  LfZ^  +  Lg 

'  L9X;  +  LjoTi  + +  1 


(39) 
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and  leave  all  the  parameters  L\...Ln  variable.  [When  interpreted  as  a  projective  transform,  the  same  11 
coefficients  are  constrained  to  only  nine  degrees  of  freedom,  corresponding  to  the  six  exterior  camera 
parameters  from  Eq.  (8)  and  the  three  interior  camera  parameters  Uc,  v^,  and  fc  from  Eq.  (7).]  Because  of 
the  two  additional  degrees  of  freedom,  a  lower  rms  value  is  obtained  in  Eq.  (38),  giving  a  (potentially) 
improved  mapping  between  the  2D  and  3D  coordinates.  In  a  more  sophisticated  approach,  additional 
terms  corresponding  to  specific  deformation  modes  of  the  test  article  are  included,  such  as  bending  and 
twisting  of  a  wing.^^’^’  However,  a  detailed  treatment  of  such  model  deformation  issues  is  beyond  the 
scope  of  the  present  study. 

4.8  IMPLEMENTATION  OF  ALGORITHM  AND  ITS  USE 

Figure  23  illustrates  a  possible  implementation  of  the  theory  from  Sections  4.3  through  4.7  to  perform 
automatic  image  registration  by  (1)  estimating  the  position  x^*  and  attitude  matrix  of  the  test  article 
with  respect  to  the  test  facility;  (2)  calculating  the  resulting  mapping  coefficients  for  transformation  of  3D 
model  coordinates  to  2D  image  coordinates;  (3)  calculating  initial  predictions  of  the  locations  of  all 
visible  registration  targets,  including  search  templates;  (4)  finding  the  exact  locations  of  the  targets;  and 
(5)  updating  the  resulting  mapping  coefficients-for  example,  by  calculating  updated  values  of  Xj*  and 
R,*. 


Tunnel  Camera 
Parameters  Parameters 


3D  Grid  of 
Test  Artiele 


Update 
Map  Coef. 

~ir- 

Done 


Figure  23.  Image  Registration  Flow  Diagram 

Four  sets  of  input  are  required:  (1)  a  list  of  model  coordinates  x„  normals  n„  and  target  diameters  c/, 
for  all  registration  targets;  (2)  exterior  parameters  Pc*  =  (xc*,  yc*,  Zc*,  (tie*,  k^*,  cOe*)  and  interior 
parameters  =  {/[,  v^,  K^}  for  each  camera;  (3)  a  set  of  fixed  tunnel-related  parameters  t  =  {ao,  x^, 

Zm,  E};  and  (4)  a  set  of  angles  &k  =  {a,  as,  it)s}  for  each  data  point.  Optionally,  a  3D  grid  of  the  test  article 
can  be  provided  as  a  fifth  input,  if  z-buffer  visibility  checks  are  to  be  performed  based  on  Eq.  (25). 
Additional  parameters  that  have  to  be  specified  are  listed  below. 

(1)  The  parameters  s,  Nu,  Ny,  and  Zw  from  Eqs.  (23)  and  (25),  which  are  used  to  verify  target  visibility. 

(2)  The  number  of  points  n  from  Eq.  (26),  used  to  calculate  the  search  templates.  (A  typical  value  is  n 
=  90,  though  smaller  values  would  probably  suffice.) 

(3)  Minimum  acceptable  correlation  values  Jmin  and  Jmin*  for  the  nonaugmented  and  augmented 
correlation  values,  respectively  (typical  values  are  J^m  =  0.5  and  Jmin*  =  0.75),  and  a  maximum  acceptable 
value  for  the  magnitudes  of  the  subpixel  shifts  Am  and  Av  from  Eq.  (35)  (typically,  0.7  pixels  each). 

(4)  The  magnitude  of  the  search  radius  around  the  initial  target  prediction  within  which  the  actual 
target  is  to  be  found  (typically,  5  to  10  pixels). 

(5)  A  minimum  acceptable  fraction  of  targets  (found  versus  predicted)  at  the  end  of  a  correlation 
search  (typically  0.9  if  a  z-buffer  is  used  to  ascertain  target  visibility). 
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(6)  A  maximum  acceptable  deviation  Ssmax  between  measured  and  fitted  image  coordinates  {ui,  v,)  for 
individual  targets  i  in  Eq.  (38)  to  help  flag  outlier  points  (typically,  Ssmax  =  2  pixels). 

(7)  A  maximum  acceptable  rms  fit  error  for  all  targets  combined  in  Eq.  (37)  (typically,  1  pixel). 

Selection  of  these  various  parameters  depends  on  the  test  being  performed  and  usually  requires  some 
empirical  optimization.  (Quoted  values  are  typical  for  PSP  testing  in  AEDC’s  16-ft  Transonic  Wind 
Tunnel  16T.) 

4.9  SAMPLE  RESULTS 

The  image  registration  algorithm  described  in  Sections  4.3  through  4.8  has  been  used  extensively  as 
the  enabling  step  for  automatic  processing  of  thousands  of  PSP  images  from  Tunnel  16T  (see,  especially. 
Ref.  24).  Performance  results  of  the  template-based  search  algorithm  are  documented  in  Ref  29.  Use  of 
the  augmented  template  technique  typically  results  in  99  percent  or  more  of  visible  registration  targets 
being  found  correctly.^^  (For  poor  lighting  conditions,  this  requires  extension  of  the  augmented  template 
technique  as  described  in  Ref  29.)  As  a  single  example.  Fig.  2  shows  a  registered  image  from  a  PSP  test 
of  an  F-16C,  with  closeups  of  some  of  the  registration  targets.  Target  diameters  in  Fig.  2  (at  about  two 
image  pixels)  were  smaller  than  ideal  (about  five  pixels),  but  they  still  allowed  target  centers  to  be 
determined  to  accuracies  of  about  0.07  image  pixels  in  both  the  u-  and  v-directions. 

The  most  critical  element  of  the  automatic  image  registration  scheme  is  the  prediction  model  from 
Section  4.4.  (Without  good  initial  predictions  for  the  target  locations,  the  simple  yet  powerful  search 
scheme  from  Section  4.6  becomes  much  harder  to  implement;  see,  e.g..  Ref  16.)  An  important  figure  of 
merit,  therefore,  is  the  worst-case  rms  prediction  error  for  visible  target  locations  in  an  image,  defined 
analogous  to  Eq.  (37),  but  with  predicted  image  coordinates  taking  the  place  of  fitted  coordinates.  Table  6 
lists  such  maximum  prediction  errors  for  three  tests  in  Tunnel  16T,  involving  an  F-16C  jet,  NASA’s  X-38 
Crew  Return  Vehicle,  and  a  Joint  Strike  Fighter  model.  The  maximum  rms  prediction  errors  are  2.4,  3.1, 
and  3.7  image  pixels,  respectively,  assuming  that  all  fixed  parameters  (such  as  camera  focal  lengths  and 
effective  sting  length)  have  been  optimized.  Average  prediction  errors  are  even  smaller,  at  around  one 
image  pixel  for  all  three  tests.  Also  shown  in  Table  6  are  the  maximum  rms  prediction  errors  that  result  if 
sting  deflection  is  not  taken  into  account,  giving  maximum  rms  errors  of  7.8,  3.8,  and  9.2  image  pixels, 
respectively.  Inclusion  of  a  droop  term  (to  account  for  bending  of  the  sting  under  the  force  of  gravity, 
typically  on  the  order  of  0.3  deg  or  less)  has  also  been  implemented,  but  this  did  not  result  in  a  significant 
reduction  of  the  maximum  prediction  errors.  Alternate  sting  positioning  models  are  also  possible.  For 
example,  some  use  a  yaw  head  or  a  knuckle-sleeve  system  that  compensates  for  vertical  motion  of  the  test 
article. 

4.10  CONCLUSION 

An  efficient  technique  for  registration  of  digital  images  of  sting-mounted  test  articles  undergoing 
rigid-body  motion  has  been  implemented  at  AEDC  specifically  to  accomplish  automatic  data  reduction  of 
PSP  images.  The  technique  relies  on  the  initial  prediction  of  the  location  of  all  registration  targets  that 
should  be  visible  in  an  image,  the  use  of  an  augmented  template  search  to  find  the  actual  target  locations 
to  subpixel  accuracy,  and  a  nonlinear,  least-squares  fitting  algorithm  to  establish  the  exact  mapping 
transformation  from  3D  model  coordinates  to  2D  image  coordinates.  In  the  process,  an  accurate 
measurement  of  the  position  and  attitude  of  the  test  article  is  obtained.  Two  elements  of  the  technique,  in 
particular,  are  key  to  its  successful  implementation  and  are  believed  to  be  novel.  These  are  the  use  of  a 
sting  deflection  model  under  air  loads  and  the  use  of  an  augmented  template  scheme  for  correlation-based 
target  finding.  Results  indicate  that  initial  prediction  of  target  locations  is  possible  to  within  about  five 
image  pixels  and  that  the  final  accuracy  of  target  localization  is  better  than  one-tenth  of  an  image  pixel. 
The  technique  has  been  used  successfully  on  thousands  of  images  acquired  in  Tunnel  1 6T  and  should  be 
applicable  to  optical  testing  in  other  test  facilities. 
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Table  6.  Selected  Data  from  Three  Wind  Tunnel  Tests 


F-16C 

X-38 

JSF 

Number  of  Cameras 

8 

8 

5 

Number  of  Images 

400 

1560 

1130 

Targets  per  Image  (Avg.) 

34 

12 

17 

Mach  Number  Range 

0  to  0.9 

Oto  1.4 

0  to  0.9 

ao,  deg 

1 

3 

8 

L,  in 

63 

88 

38 

a  Range,  deg 

0  to  26 

4  to  20 

-8  to  24 

P  Range,  deg 

0 

-4  to  4 

-8  to  8 

as  Range  deg 

- 1  to  25 

1  to  17 

-16  to  16 

(t)s  Range,  deg 

0 

-21  to  21 

-46  to  45 

Aas  Range,  deg 

2.4 

0.7 

1.5 

Max  Rms  Error  w/o  Defl,  pix 

7.8 

3.8 

9.2 

Max  Rms  Error  w/  Defl,  pix 

2.4 

3.1 

3.7 

Average  Rms  Error  w/  Defl,  pix 

0.9 

1.3 

1.2 
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5.0  EXTENSION  OF  SUBPIXEL-AUGMENTED  TEMPLATE  MATCHING  TO 
IMPERFECTLY  RENDERED  BINARY  TARGETS^ 


5.1  OVERVIEW 

In  a  previous  publication  (Ref.  1 7)  a  new  method  was  introduced  for  correlation  matching  of  binary 
targets  (e.g.,  a  light,  uniform  target  on  a  dark,  uniform  background).  This  new  method,  which  is 
summarized  in  Section  3.6,  is  based  on  the  use  of  an  augmented  search  template  that  contains  not  only  the 
occupied  pixel  fractions  of  the  synthetic  reference  that  is  being  located,  but  also  the  derivatives  of  these 
pixel  fractions  with  respect  to  shifts  of  the  template  along  the  image  axes.  In  this  section,  the  augmented 
template  technique  is  extended  to  cases  in  which  a  binary  target  is  used  for  image  registration  but  is  not 
rendered  as  a  binary  image  feature  because  of  experimental  artifacts  such  as  specular  reflection,  mixed 
image  contrast,  or  partial  occlusion.  It  is  shown  that  subpixel  localization  of  the  target  can  still  be 
accomplished,  by  extending  the  concept  of  a  single  augmented  template  to  a  composite  template,  a 
method  that  consists  of  overlapping  three-by-three  subtemplates.  Examples  of  noncomposite  and 
composite  augmented  template  matching  are  demonstrated  experimentally. 

5.2  INTRODUCTION 

Many  engineering  disciplines  require  the  (automatic)  detection  of  synthetic  references  (SRs)  in  digital 
images  for  machine  vision,  image  registration,  optical  diagnostics,  metrology,  etc.  A  common  approach  is 
the  use  of  search  templates  and  correlation  matching.  ’  ’  Typically,  the  search  template  is  stepped 
across  the  scene  in  one -pixel  increments.  Detection  of  the  SR  occurs  when  the  chosen  similarity  measure 
exceeds  a  threshold  value.  To  achieve  subpixel  resolution,  sampled  correlation  interpolation  (SCI)  may  be 
used,  ’  ■  or  multiple  search  templates  may  be  used,  each  calculated  for  a  different  subpixel  placement 
of  the  In  Ref  17  a  new  approach  is  explored,  based  on  the  insight  that  it  is  possible  to  achieve 

subpixel  resolution  by  using,  in  the  search  template,  not  only  the  occupied  pixel  fractions  (as  is  done  in 
conventional  correlation  matching),  but  also  the  derivatives  of  the  occupied  pixel  fractions  with  respect  to 
shifts  of  the  template  along  the  image  axes.  This  approach  has  two  advantages  over  conventional 
correlation  matching.  First,  it  produces  larger  correlation  values  for  unfavorable  subpixel  displacements 
between  SR  and  scene,  thereby  causing  fewer  SRs  to  go  undetected.^^  Second,  the  new  approach  allows 
the  subpixel  location  of  the  SR  to  be  determined  more  precisely  than  is  possible  with  either  SCI  or 
centroiding."^^’"^"^  In  Ref  1 7,  these  properties  of  augmented  template  matching  are  demonstrated  for  strict 
use  of  normalized  correlation  as  the  similarity  measure.  Pertinent  highlights  from  Ref  1 7  are  summarized 
in  Sections  5.3  and  5.4  (as  well  as  in  Section  3.6).  In  Sections  5.5  through  5.7  the  theory  from  Ref  17  is 
applied  in  practice  and  is  generalized  for  situations  in  which  experimental  artifacts  preclude  the  use  of 
normalized  correlation  in  its  conventional  form. 

5.3  REVIEW  OF  CONCEPT 

Let  the  term  Jfg  denote  a  suitable  similarity  measure  between  a  template  / and  a  scene  g.  For  example, 
Jfg  might  be  the  normalized  correlation  (Sections  5.4  and  5.5)  or  a  more  general  measure  (Sections  5.6  and 
5.7).  Typically,  the  template  is  calculated  for  a  fixed  relative  placement  of  the  SR,  e.g.,  with  the  SR 
centered  on  a  pixel.  A  search  is  then  performed  by  stepping  the  template  across  the  scene  in  one -pixel 
increments.  An  array  of  values  Jfg  (m,n)  is  thus  obtained,  where  m  and  n  are  integers.  Detection  of  the  SR 
occurs  when  Jfg  exceeds  a  threshold  value  J^m-  ’  ’  In  Ref  17  subpixel  resolution  is  achieved  by 


^Previously  published  as  “Extension  of  Subpixel  Augmented  Template  Matching  to  Imperfectly  Rendered  Binary 
Targets,”  in  Optical  Engineering,  Vol.  43,  No.  3,  March  2004,  pp.  639-647. 
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expanding  the  search  template  as  a  first-order  Taylor  series  in  terms  of  a  set  of  to-be-determined  subpixel 
shifts  Ax  and  Ay  according  to 

/  =  /o  +  fx^  +  fy^y-  (40) 

Figure  24  (see  also  Figs.  19  and  20)  shows  an  example  of  such  an  augmented  template  (AT)  for  a 
binary  SR,  i.e.,  an  SR  that  is  chosen  to  yield  a  high-contrast  image  feature  with  roughly  uniform 
foreground  and  background  signals.  (Contrast  with  this  the  use  of  the  term  “binary  template”  to  denote  a 
template  that  consists  of  only  ones  and  zeros."^^)  For  a  binary  SR  in  the  current  sense  (with  0  </o  ^  1),  the 
fx  and  fy  terms  can  be  obtained  as  occupied  height  and  width  differences  at  opposite  sides  of  a  pixel,  as  is 
illustrated  in  Fig.  25. 


Figure  24,  Example  of  Augmented  Template 


/  =  .43  -  0  =  .43 
/  =  .55  -0=.55 

Jy 


f^=\  -0=1 
/  =.68  -  .55  =  .13 

■'y 


Figure  25.  Determination  of  AT  Derivative  Terms 
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When  Eq.  (40)  is  substituted  into  the  chosen  similarity  measure,  the  discretely  sampled  Jfg  (m,n) 
becomes  a  function  of  the  shifts  Ax  and  Ay,  i.e.,  Jfg  =  Jfg(m+Ax,n+Ay).  The  role  of  the  AT  can  now  be 
defined  as  follows:  for  a  given  set  of  integer  shifts  (m,n),  find  the  additional  shifts  (Ax,  Ay)  that  maximize 
Jfg  (w+Ax,n+Ay).  This  maximized  correlation  value  is  denoted  with  an  asterisk,  as  shown  below. 

J jg  {m,  n)  =  max  (w  +  Ax,  n  +  Ay)].  (41) 

Ax,  Ay 

This  augmented  correlation  still  needs  to  be  optimized  as  a  function  of  the  integer  shifts  (m,n),  as  is 
the  case  in  conventional  correlation  matching.  Figure  26  gives  a  typical  implementation  of  this  concept: 
an  initial  placement  (xo,yo)  of  the  SR  is  selected,  and  the  AT  is  calculated  (Steps  1  and  2).  A  standard 
correlation  search  is  then  performed  (Steps  3  through  5),  until  Jfg  exceeds  a  preset  threshold  Jmin-  At  this 
point,  the  AT  is  invoked  (Steps  6  through  9)  to  calculate  shifts  Ax  and  Ay.  A  check  is  performed  (Step  7) 
to  verify  that  the  magnitude  of  the  shift  is  below  a  threshold  value  Asmax,  typically  one  pixel.  Definitive 
detection  of  the  SR,  at  a  location  (xo+w+Ax,  yo+n+Ay),  is  flagged  when  the  AT-based  correlation  Jfg  * 
exceeds  a  threshold  value  Jmin*  (which  is  set  higher  than  Jmm).  The  SR  being  found,  the  AT  is  recalculated 
and  recentered  iteratively  (Steps  10  through  14)  until  the  magnitudes  of  the  subpixel  shifts  (Ax,  Ay)  are 
smaller  than  their  respective  precisions.  Examples  of  this  approach  are  given  in  Sections  5.4  through  5.7, 
below. 


Figure  26.  Flow  Diagram  for  Augmented  Template  Search  and  Centering 

5,4  REVIEW  OF  THEORY  FOR  BINARY  SCENE 

When  scene  and  template  are  related  by  a  strict  gray-level  transformation,  i.Q.,  g  =  b  +  cf,  where  the 
background  b  and  contrast  c  are  constants,  a  suitable  similarity  measure  is  provided  by  the  normalized 
correlation^^’^® 

Jfg  -[i;(/-</»(g-<g»]Y[i;(/-</»2][S(g-<g»2].52^/Sy^Sgg  .  (42) 
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When  the  AT  expansion  from  Eq.  (40)  is  substituted  into  Eq.  (42),  additional  correlation  sums  arise, 
namely  Sf,,  Sjy,  5A,  S^y,  Syy,  S^g,  and  Syg  (see  Fig.  22).  For  example,  Sf,  =  Wo-{fo})(fx-{fx}),  with  {Q  =  0.  As 
is  shown  in  Ref.  17,  maximization  of  the  ensuing  Jfg  with  respect  to  the  subpixel  shifts  Ax  and  Ay  is 
facilitated  by  initial  placement  of  the  SR  such  that  the  derivative  terms  fx  and/v  are  orthogonal  to  the  pixel 
fractions  /o,  in  the  sense  that  =  Ay,  =  0.  In  this  case,  the  shifts  Ax  and  Ay  that  maximize  Jfg  are  given 
by'^: 


'Ax) 


( ^xx 

Sfg 


.  \ 

-1 

xy 

yyj 

S-'sj  . 


(43) 


The  sums  Sff  and  Sfg  in  Eq.  (43)  are  to  be  computed  using  the  pixel  fractions  fo.  The  sums  Sxg  and  Syg 
determine  the  magnitudes  of  the  subpixel  shifts  and  thus  play  a  role  somewhat  similar  to  Sobel  or  Prewitt 
edge  operators."*®  (See  also  an  older  paper  in  which  the  gradients  of  the  search  template  are  used  for 
subpixel  searches  of  image  features."*®)  Given  the  shifts  Ax  and  Ay  from  Eq.  (43),  the  augmented 
normalized  correlation  is  given  in  terms  of  the  nonaugmented  value  from  Eq.  (42)  by: 

J*j-g  -  Jfg  fg  [because  Det(S)>0].  (44) 

^fg 


In  the  presence  of  noise,  the  rms  variance  of  the  shifts  Ax  and  Ay  may  be  estimated  from  the 
computed  maximum  correlation  values  as:*^ 


5x^=  ^43 


N- 


j(l-4)(s"‘)-.  Spli 


likewise , 


(45) 


where  N  is  the  number  of  pixels  in  the  template.  Note  that  perfect  correlation  (Jy*=l)  implies  perfect 
centering.  The  results  from  Eqs.  (43)  through  (45)  can  be  generalized  for  Ay  Q,Sfy^Q  (see  Ref  17),  as  is 
the  case  during  the  iterative  centering  loop  of  Steps  10  through  14  in  Fig.  42. 

5.5  APPLICATION  TO  BINARY  SCENES 

In  this  section  the  discussion  turns  from  the  material  in  Ref.  17  to  discuss  the  application  of  AT-based 
template  matching  to  image  data  taken  in  large-scale,  transonic  wind  tunnels  at  AEDC  and  at  NASA 
Ames  Research  Center.  In  both  facilities,  subscale  models  of  flight  vehicles  are  tested  that  are  nominally 
6  ft  in  size.  Cameras  are  mounted  at  a  distance  of  about  8  ft.  Data  from  three  tests  are  considered  here  (see 
Table  7).  The  first  two  are  pressure-sensitive  paint  (PSP)  tests;  the  third  is  a  model  attitude  test. 


Table  7.  Parameters  for  Wind  Tunnel  Tests  Considered  in  This  Study 


Test  I 

Test  II 

Test  III 

Type  of  Test 

PSP 

PSP 

Model  Attitude 

Test  Article 

F-16C 

F-16C 

Model  Sting 

Test  Center 

AEDC 

NASA  Ames 

AEDC 

Imaged  SR  Definition 

Good 

Good 

Various 

Avg.  Imaged  SR  Diameter,  pixels 

2.6 

5 

13 

Number  of  SRs  in  Study 

1163 

479 

3009 

Detected  SRs  using  AT,  percent 

96.4 

92.9 

74 

Detected  SRs  using  CAT,  percent 

>99.7 

>99.7 

>99.7 
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5.5.1  Pressure-Sensitive  Paint  Tests 


Figure  27  shows  a  sample  CCD  image  from  Test  I.^  Images  from  such  a  PSP  test  are  converted  to 
surface  pressure  maps  by  exploiting  the  principle  of  oxygen  quenching  of  fluorescence  from  a  dye  in  the 
painted  surface.’  Processed  data  from  up  to  eight  cameras  are  mapped  to  a  3D  grid  of  the  test  article.  To 
accomplish  this,  small,  circular  targets  are  applied  to  the  test  article.  Automatic  image  registration  is  used 
to  track  these  targets  with  subpixel  accuracy. 


Figure  27.  Demonstration  of  AT-Based  SR  Detection  and  Centering  for  Test  I  of  Table  7 

At  the  heart  of  the  tracking  algorithm,  AT  search  templates  are  computed  for  each  SR  based  on  the 
following  five  pieces  of  information  (see  Section  4.5):  (1)  the  physical  coordinates  of  the  SR  at  a 
reference  orientation  of  the  test  article,  as  measured  using  a  coordinate-measuring  machine;  (2)  the 
physical  diameter  of  the  SR;  (3)  the  surface  normal  of  the  SR  at  the  chosen  reference  orientation;  (4)  the 
estimated  position  and  orientation  of  the  test  article  based  on  the  tunnel-supplied  pitch  and  roll  angles  of 
the  model;  and  (5)  the  calibrated  positions  and  orientations  of  the  cameras  with  respect  to  the  wind  tunnel, 
as  well  as  their  internal  parameters.  Photogrammetry  is  used  to  relate  2D  image  coordinates  to  3D 
physical  coordinates.’^’^*  Insets  in  Fig.  27  show  some  of  the  detected  SRs  for  Test  I  along  with  augmented 
correlation  values  Jfg*  and  estimated  centering  precisions  Ss  =  from  Eq.  (45),  in  image  pixels. 

In  the  wind  tunnel  application,  actual  accuracies  of  the  SR  centering  determinations  are  known  only 
as  rms  deviations  between  measured  and  fitted  SR  centers  for  all  SRs  combined  (see  Section  5.10).  These 
deviations  are  typically  in  the  range  from  0.3  to  0.8  pixels  rms  but  include  effects  of  model  deformation 
under  air  loads,  errors  in  initial  measurement  of  the  physical  coordinates  of  the  SRs,  and  camera 
calibration  errors,  if  any.  Even  so,  use  of  the  AT  technique  in  Tests  I  and  II  resulted  in  successfiil 
detection  of  96  and  93  percent,  respectively,  of  the  imaged  SRs.  (For  Test  I,  an  earlier  search  algorithm 
based  on  centroiding  achieved  only  a  92-percent  success  rate.)  Best  performance  of  the  AT  algorithm  was 
obtained,  particularly  for  small  SRs,  by  including  a  one-pixel-wide  skirt  around  the  minimum  bounding 
box  of  a  projected  SR,  as  shown  in  Fig.  24. 

Figure  28  shows  average  convergence  rates  for  the  AT-based  centering  loop  from  Fig.  26  (steps 
10  through  14)  for  Tests  I  and  II.  The  plotted  quantity  is  the  ratio  of  the  calculated  step  size,  (Ax^-l-Ay^)’^^, 
to  the  computed  precision,  (8x^-l-8y^)’^^,  at  each  iteration.  (Iteration  zero  corresponds  to  the  initial  search 
step  on  the  left  side  of  Fig.  26;  the  calculation  of  the  precision  for  this  initial  step  is  not  shown  in  Fig.  26.) 
Figure  28  shows  that  a  few  iterations  of  the  centering  loop  suffice  to  determine  the  subpixel  location  of  an 
SR  to  within  the  achievable  resolution,  given  the  noise  level  of  the  data. 


43 


AEDC-TR-06-6 


oj 

PJ 

o 


o 

(D 


(D 


Number  of  iterations 

Figure  28.  AT-Based  Convergence  Rates  of  Centering  Loop  in  Fig.  26 


5.5.2  Model  Attitude  Tests 

Figure  29  shows  an  image  taken  from  Test  III  in  which  both  black  and  retroreflecting  targets  were 
used  as  SRs.  From  the  perspective  of  image  registration,  Test  III  differs  from  Tests  I  and  II  primarily  in 
the  fact  that  imaged  SRs  suffer  from  a  range  of  deviations  from  pure  binary  imaging  because  of  specular 
reflection,  mixed  image  contrast,  and  occlusions.  Using  the  AT  scheme,  only  74  percent  of  SRs  were 
detected.  Sections  5.6  and  5.7  demonstrate  how  the  AT  technique  can  be  extended  to  handle  such  cases 
more  effectively. 


Figure  29.  Image  from  Test  III  of  Table  7,  with  SRs  Detected  by  the  AT  Scheme 

Figure  30  shows  an  image  of  an  18-in.  calibration  plate  that  was  used  as  part  of  Test  III.  For  this 
plate,  physical  coordinates  of  the  SRs  were  known  with  high  accuracy  (better  than  0.1  mm),  and  model 
deformation  is  not  an  issue.  In  this  case,  rms  deviations  between  measured  (i.e.,  AT-based)  image 
coordinates  and  fitted  coordinates  are  about  0.05  image  pixels  (see  Table  8).  Table  8  also  shows  good 
agreement  between  the  predicted  precisions  from  Eq.  (45)  and  the  rms  fit  errors  from  Section  5.10.  This 
result  shows  (as  is  also  seen  in  Section  5.7)  that  Eq.  (45)  is,  indeed,  a  useful  measure  of  the  accuracy  of 
SR  image  coordinates  at  the  end  of  the  centering  loop  in  Fig.  26. 

5.6  EXTENSION  TO  NONBINARY  SCENES:  THEORY 

To  apply  the  AT  technique  to  imperfectly  rendered  binary  SRs,  allowance  must  be  made  for  the  fact 
that  there  is  no  single  gray-level  transformation  for  the  SR  as  a  whole."^’  Flowever,  the  template  may  be 
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5s  =  0.04 

Figure  30.  Image  of  Calibration  Plate  in  Test  III  with  SRs  Detected  by  the  AT  Scheme 


Table  8.  Comparison  of  Centering  Errors  for  Calibration  Plate  in  Test  III  of  Table  7 


Camera  1 

Camera  2 

Avg.  Image  Contrast  (%  of  Peak  Signal) 

High  (76%) 

High  (95%) 

Avg.  Imaged  SR  Diameter,  pixels 

6 

9 

Avg.  Predicted  Precision  (x,  y  combined),  pixels 

0.06 

0.06 

Avg.  Rms  Fit  Error  (x,  y  combined),  pixels 

0.05 

0.04 

divided  into  subtemplates  (see  Fig.  31),  each  of  which  satisfies  (at  least  approximately)  a  gray-level 
transformation.  Now  a  composite  sum  of  squared  differences  (CSSD,  cf  Ref  48)  between  template  and 
scene  can  be  defined  as  follows  (see  Section  5.1 1): 

CSSD  =  2*(%^te-4j- 

where  the  sum  over  k  extends  over  each  of  the  subtemplates  that  make  up  the  template  as  a  whole.  The 
correlation  sums  Sk/f,  St/g,  and  Stgg  are  defined  as  in  Eq.  (42),  but  for  each  subtemplate  k  individually.  For 
example,  Skgg  =  l.(g-(g)kY,  where  averaging  and  summation  extend  only  over  that  part  of  the  scene  that 
corresponds  to  a  particular  subtemplate.  From  Eq.  (46),  the  composite  normalized  correlation  can  be 
computed  as 

fg  =  ^  “  CSSD  /  'Y.k  ^kffSkgg  =  Za;  ^kfg  l^k  ^¥f^kgg  >  (47) 

with  0  <  J/g  <  1.  It  is  shown  in  Section  5.11  that  the  CSSD  from  Eq.  (46)  (and  hence  the  ensuing 
correlation)  is  especially  sensitive  to  pixels  near  the  edge  of  the  SR. 
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Figure  31.  Division  of  Search  Template  into  Sub  templates 
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To  apply  the  AT  concept  to  a  composite  template,  it  is  possible  to  follow  the  approach  taken  in 
Section  5.4  and  substitute  the  AT  expansion  from  Eq.  (40)  directly  into  Eq.  (47).  However,  a  more 
tractable  problem  is  obtained  by  substituting  the  AT  expansion  into  Eq.  (46)  and  minimizing  the  CSSD 
instead.  This  leads  to  subpixel  shifts  Ax  and  Ay  given  by 


'Ax) 

( D  D  ^ 

^xx  ^xy 

-1 

[Dfx] 

Dyy  y 

D 


/ 


(48) 


where 

^  na  —  ~  ^kr)p^kap)i  P  —  fr,X,y,  q  —  X,y.  (49) 

The  subtemplate  sums  S/ag,  S^g,  Skxx,  S/ay,  Skyy,  Skjk,  and  Skfy  are  defined  analogously  to  their 
counterparts  in  Section  5.4,  e.g.,  Skfx  =  '^(fo-{fo)k)(fx-{fx)k),  with  {Qk  0  in  this  case.  The  augmented 
correlation  value  is  now  given  in  terms  of  the  nonaugmented  value  from  Eq.  (47)  by  [cf.  Eq.  (44)]: 

^fg  =  ^fg  +  (‘I/  •  )/X^ SkffSkgg  >  Jfg  [because  Det(D)  >  0] .  (50) 

In  the  presence  of  noise,  the  rms  variance  of  the  shifts  Ax  and  Ay  may  again  be  estimated  from  the 
augmented  correlation  values  [cf.  Eq.  (45),  though  additional  approximations  are  required  in  the  present 
case],  yielding: 

g^2  ^  (l  -  )(d-i  ,  5y2  likewise ,  (51) 

with  N  being  the  number  of  pixels  in  the  template.  To  distinguish  results  based  on  the  use  of  a  composite 
template  from  those  based  on  a  noncomposite  template,  the  acronym  CAT  (composite  AT)  is  introduced 
for  the  former. 

5.7  EXTENSION  TO  NONBINARY  SCENES:  PRACTICE 

Different  tiling  strategies  were  investigated  for  subdividing  search  templates  in  cases  where 
correlation  matching  based  on  the  AT  scheme  from  Sections  5.4  and  5.5  failed.  Of  these  tilings,  frilly 
overlapping  3  by  3  subtemplates  (the  case  shown  in  Fig.  31)  were  found  to  give  a  satisfactory 
compromise  between  ease  of  implementation  and  performance,  though  this  tiling  choice  may  not  be 
optimal  for  all  applications.  Redundancy  introduced  by  selecting  overlapping  subtemplates  appeared  to 
benefit  detection  of  poorly  defined  SRs,  though  at  the  expense  of  higher  computational  cost  as  compared 
to  nonoverlapping  3  by  3  subtemplates.  Subtemplates  of  size  2  by  2  were  found  to  give  unsatisfactory 
performance  for  poorly  rendered  SRs. 

Implementation  of  the  CAT  scheme  was  accomplished  by  a  small  modification  of  the  flowchart  from 
Fig.  26.  Initially,  the  SR  search  is  performed  using  the  noncomposite  AT  scheme  from  Sections  5.4  and 

5.5.  Only  if  this  scheme  fails  to  detect  the  desired  SR  is  the  CAT  scheme  applied.  This  is  done  by 
repeating  Steps  3  through  9  in  Fig.  26,  using  the  CAT-based  correlation  expressions  for  J/g  from  Section 

5.6.  Once  the  SR  is  found,  the  centering  loop  in  Steps  10  through  14  of  Fig.  26  is  performed  for  the 
scheme  that  successfully  produced  the  SR,  i.e.,  either  AT  or  CAT.  Typical  values  for  the  threshold 
parameters  from  Fig.  26  were  Jmin  =  0.50,  Asmax  =  1  (pixel),  and  Jmin*  =  0.75. 

Applying  the  CAT  scheme  to  the  data  from  Table  7  increased  the  fraction  of  detected  SRs  to  above 

99.7  percent  for  each  of  Tests  I,  II,  and  III.  Especially  for  Test  III,  this  is  a  significant  improvement  over 
the  AT  scheme  (which  detected  only  74  percent  of  all  targets,  as  shown  in  Table  7).  Figure  32  shows 
targets  from  Test  111  that  were  detected  using  the  CAT  scheme  (but  missed  by  the  AT  scheme),  with 
associated  correlation  values  and  estimated  localization  precisions.  Even  an  SR  that  is  occluded  almost 
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completely  (Fig.  32d)  is  found  correctly  by  the  CAT  algorithm.  Likewise,  Fig.  33  shows  a  few  SRs  from 
Tests  I  and  II  that  were  detected  by  CAT  but  missed  by  AT  because  of  extremely  oblique  viewing  angles. 
Convergence  rates  for  CAT  centering  was  found  to  be  similar  to  the  AT  case  in  Fig.  28. 


4*=0.83,  8s=0.07 
a. 


J,g*=0.77,  8s=0.07  8.v=0.06  J,  *=0.1A,  8.v=0. 1 0 

b.  c.  d. 


4*=0.75,  8.s=0.08 
e. 


Figure  32.  Detected  SRs  from  Test  III  Found  Using  the  CAT  Scheme 


4*=0.94,  8.v=0.06  4*=0.94,  8s=0. 1 0 

a.  Test  1  b.  Test  11 


Figure  33.  Detected  SRs  from  Tests  I  and  II  Found  Using  the  CAT  Scheme 

Figures  34  through  36  show  further  comparisons  between  AT  and  CAT  for  SRs  that  were  detected 
using  both  schemes.  Data  points  represent  average  values  taken  over  all  detected  SRs.  As  expected  (Fig. 
34),  CAT  generally  produces  higher  correlation  values  than  does  AT,  particularly  for  imperfectly 
rendered  SRs  (for  example,  SRs  with  AT  correlations  below  about  0.95).  Still  (see  Fig.  35),  the  calculated 
precisions  of  the  resulting  centering  determinations  based  on  Eqs.  (45)  and  (51)  are  roughly  equal  for  the 
two  schemes,  indicating  that  the  calculated  precisions  are  mutually  consistent. 


0.75  0.80  0.85  0.90  0.95  1.00 

Jf*  (AT) 

Figure  34.  Comparison  of  AT  and  CAT  Correlation  Values 
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5^  (AT),  pixels 

Figure  35.  Comparison  of  Predicted  Centering  Precisions  for  Data  from  Fig.  34 


^  0.00  0.05  0.10  0.15  0.20 

(55^.j.^+  pixels 

Figure  36.  Comparison  of  Measured  and  Predicted  Precisions  for  Data  from  Fig.  34 

Actual  (“truth”)  image  centers  were  not  available.  Nevertheless,  to  assess  how  well  Eqs.  (45)  and  (51) 
predict  the  actual  centering  errors,  the  magnitude  of  the  difference  between  AT-based  and  CAT-based 
image  centers  is  plotted  in  Fig.  36  as  a  function  of  the  estimated  rms  difference  based  on  those  two 
equations.  A  theoretical  simulation  with  Gaussian  noise  (the  dashed  line  marked  “Theory”  in  Fig.  36) 
suggests  a  slope  of  0.885  between  these  two  quantities.  Observed  slopes  for  Tests  I,  II,  and  III  are  roughly 
0.5,  1.4,  and  2.1,  respectively.  Figure  36  thus  suggests  that  the  predicted  precisions  based  on  Eqs.  (45) 
and  (51)  provide  a  reasonable  estimate  of  the  actual  accuracy  of  the  calculated  target  centers. 

5.8  CONCLUSION 

The  augmented  template  scheme,  developed  in  Ref  1 7  for  binary  image  features,  has  been  demonstrated 
successfully  in  practice,  showing  that  it  is  possible  both  to  find  the  centers  of  SRs  to  subpixel  accuracy 
and  to  estimate  the  precisions  of  the  computed  center  coordinates  from  the  associated  correlation  values. 
In  particular,  it  has  been  demonstrated  that  the  AT  technique  can  be  extended  to  imperfectly  rendered 
binary  features  by  dividing  the  search  template  into  a  set  of  overlapping  subtemplates,  leading  to  a 
composite  AT  (i.e.,  a  CAT)  scheme.  Using  a  combination  of  AT  and  CAT  search  templates,  a  detection 
rate  of  better  than  99.7  percent  was  demonstrated  for  circular  targets  that  are  used  for  automatic  image 
registration  for  pressure-sensitive  paint  measurements  and  model  attitude  measurements  on  scale  models 
of  flight  vehicles  in  wind  tunnels.  This  detection  rate  is  particularly  impressive  for  those  cases  in  which 
experimental  artifacts  such  as  specular  reflection,  mixed  image  contrast,  and  occlusions  preclude  the  use 
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of  normalized  correlation  per  se.  The  AT  and  CAT  schemes  should  be  broadly  applicable  in  other 
applications  in  which  synthetic  references  are  used  for  optical  tracking. 
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5.10  APPENDIX  A:  GLOBAL  RMS  FIT  ERROR 

Let  (x„  y,)  denote  one  of  N  visible  target  locations  in  the  image.  To  check  whether  the  N  image 
coordinates  (x„  y,)  are  consistent  with  the  corresponding  rigid-body  coordinates  of  the  targets  (which  are 
measured  initially  using  a  coordinate-measuring  machine),  a  projective  transform  fit  of  the  N  target 
locations  is  performed.’^  Let  (m„  v,)  denote  the  fitted  image  coordinates.  The  rms  fit  error  can  now  be 
defined  as 


r  N 

7  7ll 

{Xi  -UiY  +{yi  -V;)  Jj 

(52) 


This  rms  fit  error  (expressed  in  image  pixels)  serves  as  a  figure  of  merit  for  all  targets  combined, 
though  it  includes  the  effects  of  imperfect  pretest  determination  of  the  3D  rigid-body  coordinates  of  the 
targets,  deformation  of  the  test  article  caused  by  air  loads,  and  camera  calibration  errors,  if  any. 


5,11  APPENDIX  B:  JUSTIFICATION  OF  COMPOSITE  SIMILARITY  MEASURE 

To  justify  the  use  of  Eq.  (46)  as  a  similarity  measure,  one  can  divide  the  template  into  a  set  of 
subtemplates  k  and  assume  that  each  is  rendered  by  a  strict  gray-level  transformation  into  the  scene, 
according  to  g-{g)k  =  Ck  (f-{f)k),  where  ct  is  the  local  contrast  of  the  k-th  subtemplate.  The  following  (with 
Wk  to-be-defmed  weight  factors  and  pixel  sums  ij  extending  over  subtemplates  k  individually)  is  then  a 
measure  of  the  mismatch  between  scene  and  template: 

CSSD  -  X Z  [gij  -  {gij  )k-ck{fij-  {fij)k)} 

k  ijGk 

=  Z  ^^kgg  ~  2cA; ^kfg  +  ^kff  ]  • 
k 


For  fixed  weights  Wk,  the  CSSD  from  Eq.  (53)  is  minimized  when  q  =  SkfJSkff-,  leading  to 

CSSD  =  2;|i-[S  ]. 

k 


(54) 


To  obtain  an  expression  that,  when  substituting  the  AT  expansion  from  Eq.  (40),  gives  rise  to  a  linear 
minimization  problem  in  the  subpixel  shifts  Ax  and  Ay,  it  is  desirable  to  remove  the  factors  from  the 
denominators  in  Eq.  (54).  This  is  accomplished  by  setting  the  weights  Wk  equal  to  Skfj,  leading  to  the 
CSSD  from  Eq.  (46).  This  choice  of  weights  has  the  additional  advantage  that  the  composite  similarity 
measure  from  Eqs.  (53)  and  (46)  is  weighted  preferentially  toward  parts  of  the  SR  with  large  values  of 
Skff,  (i.e.,  toward  subtemplates  near  the  edge  of  the  SR). 
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6.0  IMPROVED  DATA  PROCESSING  FOR  PRESSURE-SENSITIVE  PAINT 
MEASUREMENTS  IN  AN  INDUSTRIAL  FACILITY  * 


6.1  OVERVIEW 

In  this  section,  six  improvements,  which  give  rise  to  eleven  benefits,  are  described  that  have  been 
implemented  at  AEDC  to  improve  the  performance  of  the  processing  scheme  from  Section  2.0  of  this 
report  for  data  reduction  of  PSP  images  from  AEDC’s  16-ft  Transonic  Wind  Tunnel  (16T).  The 
improvements  involve  a  modification  of  the  order  of  the  main  processing  steps,  improved  mapping  of  2D 
image  data  to  the  3D  grid  of  the  test  article,  averaging  camera  views  as  opposed  to  selecting  the  best 
view,  and  the  implementation  of  two  additional  diagnostics  (pixel-based  and  camera-based)  to  assess  the 
precision  of  the  processed  data.  Results  of  the  improvements  are  illustrated  for  a  PSP  test  on  a  scale 
model  of  the  NASA  Space  Shuttle  in  launch  configuration,  in  which  AEDC’s  new,  eight-camera, 
lifetime-based  PSP  system  was  used  to  acquire  the  data. 

6.2  INTRODUCTION 

Pressure-sensitive  paint  (PSP)  is  a  technology  that  is  finding  increased  use  for  measuring  full-field 
pressure  distributions  on  test  articles  in  wind  tunnels.  In  short,  a  luminescent  paint  is  applied  to  the  test 
article,  which  is  illuminated  with  short-wavelength  (for  example,  blue)  light  and  fluoresces  at  a  longer 
wavelength  (for  example,  red).  Because  of  oxygen  quenching  of  excited  probe  molecules  in  the  paint,  the 
luminescence  signal  is  dependent  on  pressure,  thereby  forming  the  basis  for  a  pressure  measurement. 

The  principles  of  PSP  have  been  reviewed  extensively  elsewhere  (see,  for  example.  Ref  1),  and  the 
first  book  dedicated  to  the  subject  has  recently  appeared."^^  This  is  not  to  say  that  PSP  has  reached  the  end 
of  its  development  cycle.  Indeed,  efforts  are  still  ongoing  around  the  world  (e.g.,  in  Canada,^®  Japan,^'’^^ 
Korea,^^  United  Kindom,^"^  Germany,^^  France,^®  as  well  as  in  the  U.S.^^'®’)  to  improve  the  technique  and 
expand  its  areas  of  applicability. 

Among  these  efforts  is  the  one  described  in  this  paper:  to  further  improve  data  processing  of  PSP 
images  from  a  large-scale  facility  such  as  Tunnel  16T  at  AEDC.  Much  progress  has  already  been  reported 
in  this  regard,’*’®^  but  it  is  shown  here  that  additional  improvements  have  still  been  possible. 

The  organization  of  this  section  is  as  follows:  In  Section  6.3,  the  status  of  PSP  data  processing  prior 
to  the  current  effort  is  reviewed.  Section  6.4  discusses  six  improvements  to  this  prior  status,  giving  a  total 
of  eleven  new  benefits,  including  two  new  diagnostic  measures  to  judge  the  precision  of  a  PSP 
measurement.  Section  6.5  illustrates  the  results  of  these  improvements  for  a  PSP  test  on  a  scale  model  of 
the  NASA  Space  Shuttle.  Section  6.6  concludes  with  a  discussion. 

6.3  REVIEW  OF  PRIOR  PROCESSING  SCHEME 

In  essence,  a  PSP  measurement  consists  of  a  measurement  of  two  or  more  signals  whose  ratio  or 
combination  of  ratios  defines  the  pressure  by  way  of  an  appropriate  calibration  function.  In  the  case  of  the 
lifetime-based  PSP  system  employed  at  AEDC,  the  two  signals  are  both  acquired  following  pulsed 
illumination  of  the  paint,  but  they  are  acquired  at  different  relative  delays  with  respect  to  the  illumination 
pulses.*’^  Eet  the  two  signals  be  denoted  Si  and  S2,  with  ratio  [See  also  Eq.  (1)]: 


R  =  S,/S^  . 


(55) 


**W.  Ruyten  and  M.  Sellers,  “Improved  Data  Processing  for  Pressure-Sensitive  Paint  Measurements  in  an 
Industrial  Facility,”  AIAA  Paper  2006-1042,  44*  AIAA  Aerospace  Sciences  Meeting  and  Exhibit,  Reno,  Nevada, 
January  9-12,  2006. 
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The  pressure,  P,  is  then  obtained  from  this  ratio  by  means  of  a  ealibration  function, /(i?),  so  that 

P  =  f{R)  =  ^hR'.  (56) 

k=0 

Typically,  as  indicated,  the  calibration  function  /  is  a  polynomial  whose  coefficients,  ft,  are 
determined  by  either  a  priori  or  in  situ  calibration  (or  a  combination  of  the  two).  Typical  numbers  of 
coefficients,  K,  are  three  to  five.  The  coefficients  fk,  in  turn,  depend  on  temperature  and  are  computed  on 
the  basis  of  the  best  available  estimate  of  the  temperature  of  the  surface  of  the  test  article. 

Calculation  of  the  signal  ratios  and  conversion  of  the  ratios  to  pressure  by  way  of  a  polynomial 
calibration  equation  are  mathematically  simple,  even  given  the  fact  that  these  calculations  must  be 
performed  on  images.  Much  harder  is  the  challenge  of  automating  the  process  whereby  image  data  are 
mapped  to  a  three-dimensional  (3D)  grid  of  the  test  article.  However,  significant  progress  has  been  made 
in  this  area,  allowing  online  processing  of  PSP  image  data  from  up  to  eight  cameras  to  be  performed 
within  ten  seconds  of  acquiring  the  data.’^’^^ 

Figure  37  illustrates  key  steps  in  the  overall  data  reduction  process:  1)  Performing  a  smoothing 
operation  on  the  raw  image  data;  2)  taking  the  ratio  of  two  images  for  each  camera;  3)  converting  the  ratio 
to  a  pressure  by  application  of  the  calibration  function  from  Eq.  (56);  4)  mapping  the  resulting  pressure 
distribution  to  the  3D  grid  for  each  camera  individually;  and  5)  merging  data  from  multiple  cameras  to 
obtain  the  final  result.  In  this  merge  process,  the  best  (most  normal)  camera  view  is  selected  for  each 
point  on  the  surface,  and  calibration  offsets  are  applied  between  cameras,  to  achieve  optimum  matching 
of  calculated  pressures  along  boundary  lines  between  adjacent  camera  views. 
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Figure  37.  Key  Processing  Steps  Prior  to  Improvements 

Not  shown  in  Fig.  37  are  1)  subtraction  of  black  images  to  adjust  for  nonzero  signals  at  zero 
illumination;  2)  thresholding  of  the  signals  before  calculation  of  the  image  ratios;  3)  image  registration, 
whereby  the  mapping  transformation  is  established  between  2D  image  coordinates  and  3D  model 
coordinates;  4)  interpolation  of  pressure  data  across  areas  of  missing  data  (particularly  across  registration 
targets);  5)  in  situ  correction  of  the  measured  PSP  pressures  based  on  measurements  from  discrete  sensors 
embedded  in  the  test  article;  and  6)  self-illumination  correction  for  reflected  emission  light.^* 

The  reason  these  additional  steps  are  not  shown  in  Fig.  37  is  not  that  they  are  not  important,  but  rather 
to  focus  the  present  discussion  on  those  steps  in  the  processing  sequence  to  which  improvements  were 
made,  as  discussed  in  Section  6.4. 

Also  indicated  in  Fig.  37  is  the  calculation  of  a  diagnostic  measure,  “rmstap,”  defined  as  the  rms 
difference  between  the  PSP  image  data  in  the  vicinity  of  pressure  taps  and  the  pressures  measured  by  the 
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taps.  So  far,  this  has  been  one  of  the  few  diagnostics  available  to  judge  the  accuracy  of  a  PSP 
measurement.  Two  additional  diagnostics  are  introduced  in  Section  6.4. 

6.4  IMPROVEMENTS  TO  PROCESSING  SCHEME 

Figure  38  shows  the  new  processing  sequence,  showing,  as  in  Fig.  37,  only  those  steps  to  which 
changes  were  made  relative  to  the  scheme  discussed  in  Fig.  37.  The  new  sequence  of  steps  is  as  follows: 
1)  Signal  ratios  are  calculated  prior  to  smoothing;  2)  the  ratios  are  mapped  to  the  3D  grid  prior  to 
conversion  to  pressure,  while  smoothing  the  data  by  interpolation;  3)  ratios  from  all  cameras  are  merged 
onto  a  single  3D  grid;  and  4)  the  ratios  on  the  3D  grid  are  converted  to  pressures  based  on  the  calibration 
function  from  Eq.  (56).  In  addition,  two  additional  diagnostics  are  computed  in  Fig.  38,  and  are  denoted 
“rmspix”  and  “rmscam.” 
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Figure  38.  Improved  Processing  Sequence 

The  changes  in  processing  between  Figs.  37  and  38  constitute  six  improvements  and  provide  a  total  of 
eleven  benefits.  These  are  explained  in  detail  in  the  following  paragraphs 

Improvement  #1:  Smoothing  After  Ratioing 

Typical  smoothing  operations  in  the  past  have  involved  taking  a  3  by  3  or  5  by  5  pixel  running 
average  (with  outlier  rejection)  on  the  black-subtracted  images  to  reduce  the  effects  of  photon  shot  noise. 
This  process  was  satisfactory  for  the  old  intensity-based  PSP  system,^  but  it  was  found  to  introduce 
additional  noise  for  the  current  lifetime-based  system.  The  reason  for  this  seems  to  be  that  the  current 
cameras  (see  Section  6.5)  have  significant  pixel-to-pixel  gain  variations.  Improvement  #I  accrues  the 
following  benefit: 

Benefit  #1:  Ratioing  before  smoothing  compensates  more  fully  for  pixel-to-pixel  gain  variations  of 
the  cameras,  and  thus  reduces  the  noise  level  in  the  processed  data  (by  about  a  factor  of  three). 

Improvement  #2:  Smoothing  as  Part  of  Mapping 

Rather  than  performing  the  smoothing  operation  immediately  following  calculation  of  the  signal 
ratios,  smoothing  is  now  implemented  in  the  form  of  an  interpolation  scheme  that  is  an  integral  part  of  the 
mapping  process.  This  scheme  is  explained  in  Fig.  39,  which  shows  a  small  part  of  the  3D  grid  projected 
(for  a  given  camera  view)  onto  the  2D  array  of  image  pixels.  To  obtain  the  mapped  image  value  (in  the 
case  of  Fig.  38,  a  signal  ratio)  at  a  grid  point  g,  the  values  from  pixels  surrounding  the  projected  grid  point 
g  are  interpolated  onto  the  grid  point,  as  explained  in  Section  6.  This  improvement  accrues  the  following 
benefits: 
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Benefit  #2:  Smoothing  at  the  mapping  stage  is  faster  than  smoothing  image  data  prior  to  mapping 
beeause  there  are  typically  fewer  points  on  the  3D  grid  than  there  are  pixels  in  the  images. 

Benefit  #3:  Smoothing  at  the  mapping  stage  in  the  form  of  interpolation  preserves  subpixel  resolution 
of  the  mapped  data,  as  explained  in  Section  6.8. 

Benefit  #4:  Smoothing  at  the  mapping  stage  allows  the  use  of  depth-buffer  masking  to  exclude  image 
pixels  that  are  close  to  a  projected  grid  point  in  2D  image  space,  but  not  in  3D  model  space.  (A  depth 
buffer  is  calculated  as  an  integral  part  of  the  mapping  transformation  even  in  the  old  processing  scheme.) 
In  particular,  depth-buffer  masking  prevents  interpolation  across  the  edges  of  the  test  article  in  the  image. 

Benefit  #5:  While  not  yet  implemented,  smoothing  at  the  mapping  stage  should  allow  for  optimum 
matching  of  the  projected  grid  resolution  to  the  pixel  spacing  in  the  image.  That  is,  a  coarse  grid  (in 
projection)  would  justify  a  larger  number  of  pixels  in  the  smoothing  operation  than  would  a  fine  grid. 


Figure  39.  Schematic  of  Mapped  Grid  Projected  onto  Image 


Improvement  #3:  Calibration  After  Mapping  During  Merge 

In  Fig.  38,  the  conversion  of  signal  ratios  to  pressures  by  way  of  the  calibration  function  from 
Eq.  (56)  is  not  performed  until  after  the  data  (in  this  case,  signal  ratios)  have  been  mapped  to  the  3D  grid 
for  each  camera,  and  ratios  from  multiple  cameras  have  been  merged  onto  a  single  3D  grid.  This 
improvement  accrues  the  following  benefits: 

Benefit  #6:  The  problem  of  performing  in  situ  calibrations  for  camera  views  that  involve  few  or  no 
pressure  taps  is  avoided. 

Benefit  #7:  A  single  calibration  can  be  applied  for  the  entire  surface  of  the  test  article.  In  particular, 
there  is  no  need  to  calculate  calibration  offsets  for  overlapping  cameras,  as  is  indicated  in  Fig.  37  by  the 
term  “match.” 

Improvement  #4:  Average  vs.  Best  View  in  Merging 

When  parts  of  the  test  article  are  imaged  by  more  than  one  camera,  the  question  arises  of  how  data 
from  multiple  cameras  need  to  be  handled.  Traditionally,  data  from  only  a  single  view  has  been  retained, 
namely  the  one  most  normal  to  the  surface.  In  the  case  of  the  AEDC  PSP  system,  this  is  a  significant 
issue  because  eight  cameras  are  typically  employed  in  the  test  section,  as  shown  in  Fig.  40.  This  means 
that  many  parts  of  the  test  article  are  viewed  by  three  or  four  cameras. 
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Figure  40.  Schematic  of  Eight-Camera  PSP  System 

In  the  new  scheme,  data  from  overlapping  cameras  are  averaged  on  the  basis  of  viewing  angle,  as 
explained  in  Section  6.9.  From  this  improvement  accrue  the  following  benefits: 

Benefit  #8:  The  effects  of  photon  shot  noise  are  reduced  because,  effectively,  more  photons  are 
collected  in  the  measurement  scheme. 

Benefit  #9:  Systematic  errors  in  camera  overlap  regions  (though  small — see  Section  IV)  are  smoothed 
out,  practically  removing  the  appearance  of  seams  from  the  mapped  data. 

Improvement  #5;  Pixel-Based  Measurement  Precision 

As  mentioned  at  the  end  of  Section  6.3,  the  sole  diagnostic  that  has  typically  been  reported  for  PSP 
data  in  the  past  is  the  rms  difference  between  the  PSP  image  data  in  the  vicinity  of  pressure  taps  and  the 
pressures  measured  by  the  taps.  In  Fig.  38,  two  additional  diagnostics  are  indicated.  The  first  of  these  is  a 
measure  of  the  effect  of  pixel-to-pixel  noise  in  the  images.  Section  6.8  explains  how  this  pixel-based 
diagnostic,  rmspix,  is  computed.  This  additional  diagnostic  accrues  the  following  benefit: 

Benefit  #10:  The  rmspix  diagnostic  provides  a  measure  of  the  precision  of  the  PSP  data.  In 
particular,  it  reveals  how  closely  the  noise  level  approaches  the  theoretical  limit,  as  determined  by  photon 
shot  noise.  This  is  discussed  further  in  Section  6.5. 

Improvement  #6:  Camera-Based  Measurement  Precision 

The  second  additional  diagnostic  in  Fig.  38  involves  the  calculation  of  the  level  of  disagreement 
between  multiple  camera  views  for  those  points  on  the  surface  of  the  test  article  that  are  viewed  by  more 
than  one  camera.  Section  6.9  explains  how  this  camera-based  diagnostic,  rms  cam,  is  computed.  This 
additional  diagnostic  accrues  the  following  benefit: 

Benefit  #11:  The  rms  cam  diagnostic  provides  a  check  on  the  implicit  assumption  in  PSP 
measurements  that  the  physical  process  by  which  excitation  light  is  converted  to  emission  light  does  not 
have  a  view-angle  dependence,  and  that  no  other  view-dependent  effects  (such  as  a  combination  of 
specular  reflection  of  excitation  light  and  spectral  leakage  in  the  measurement  process^"^)  remain 
unaccounted  for. 

Before  moving  on  to  a  discussion  of  experimental  results  in  Section  6.5,  it  should  be  noted  that  all  of 
the  improvements  discussed  in  Section  6.4  might  be  applied,  in  principle,  to  an  intensity-based  PSP 
measurement.  However,  one  of  the  practical  limitations  of  an  intensity-based  system  is  that  the  test  article 
experiences  a  model  shift  between  wind-off  and  wind-on  measurements  because  of  sting  bending  under 
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aerodynamic  loads. This  limitation  necessitates  a  separate  alignment  step  between  the  two  images  that 
make  up  the  ratio  in  Eq.  (55).  Even  then,  residual  errors  remain  because  a  model  shift  produces  changes 
in  the  illumination  distribution.  By  contrast,  in  the  current  lifetime-based  PSP  system,  the  image  pairs  that 
make  up  the  signal  ratios  in  Eq.  (55)  are  both  obtained  at  run  conditions,  and  the  complications  caused  by 
model  shifts  are  avoided. 

A  final  comment  pertaining  to  the  changes  in  processing  scheme  concerns  the  software  that  was  used 
to  implement  the  changes.  Prior  to  this  effort,  most  PSP  data  processing  at  AEDC  was  performed  by 
Green  Boot,  a  program  developed  jointly  by  NASA  Ames  and  McDonnell-Douglas  Aerospace  (i.e., 
before  MDA’s  acquisition  by  Boeing).  Currently,  however,  most  processing  of  PSP  data  at  AEDC  is 
performed  by  the  AEDC-developed  program  BlackBox,  so  named  for  its  ability  to  process  data  with 
minimal  (if  any)  user  interaction,  after  some  initial  calibrations  have  been  performed.  All  of  the  changes 
described  in  Section  6.4  were  made  in  BlackBox.  However,  the  changes  described  are  not  really  software 
specific  and  should  be  amenable  to  other  software  packages. 

6.5  SAMPLE  RESULTS  FOR  NASA  SPACE  SHUTTLE  TEST 

Most  of  the  improvements  in  processing  methodology  described  in  Section  6.4  were  implemented  to 
support  a  lifetime-based  PSP  test  of  a  three-percent  scale  model  of  the  NASA  Space  Shuttle  in  launch 
configuration  (i.e.,  with  external  tank  and  solid  rocket  boosters).  This  test  was  performed  in  October  2004 
in  the  16-ft  Transonic  Wind  Tunnel  at  AEDC.  Results  of  this  test  have  been  reported  elsewhere.®^  Here 
we  report  only  on  those  elements  of  the  test  that  highlight  the  improvements  in  processing  methodology. 

Image  data  were  acquired  with  eight  Photometries  CoolSNAP  K4  interline  transfer  cameras,  each 
with  a  spatial  resolution  of  2048  by  2048  image  pixels  and  a  gray-scale  resolution  of  12  bits,  i.e.,  0  to 
4095  counts  per  pixel.  (For  an  interesting  discussion  of  camera  selection  for  a  PSP  measurement  system, 
see  Ref  67.)  Processed  data  were  mapped  to  a  3D  grid  consisting  of  332  zones  and  643,555  grid  points. 
More  than  1,100  pressure  taps  were  available  to  compare  results  from  the  PSP  analysis  to  conventional 
measurements.  PtTFPP  in  FIB  was  used  as  the  paint.^  Pulsed  illumination  was  provided  by  banks  of 
EEDs  emitting  at  464  nm,  pulsing  every  600  ps,  and  luminescence  was  filtered  at  650  nm.  Images  were 
acquired  by  integrating  up  to  30,000  pulses  at  each  gate.  About  70  circular  registration  targets  were 
applied  to  the  test  article  to  perform  image  registration.  Data  were  processed  on  a  Einux  cluster  consisting 
of  a  master  and  eight  slaves,  each  with  two  processors.  Processing  times  were  on  the  order  of  six  seconds 
per  test  condition.  Further  details  on  the  AEDC  lifetime-based  PSP  system  are  discussed  in  Ref  68. 

Figure  41  shows  an  example  of  a  processed  data  point  for  the  test  condition  indicated.  The  display 
color  represents  pressure,  ranging  from  low  (blue)  to  high  (red). 

Figure  42  shows  the  corresponding  pixel-based  precision  rmspix,  calculated  as  shown  in  Section 
44.  Blue  indicates  small  values  of  this  rms  error,  while  red  indicates  high  values.  The  distributions  in 
Figs.  41  and  42  are  clearly  correlated,  with  better  precisions  corresponding  to  low-pressure  areas  on  the 
model  (i.e.,  areas  where  the  emission  signal  is  the  highest).  This  is  the  expected  behavior  for  a 
measurement  system  in  which  photon  shot  noise  is  the  limiting  factor. 

Figures  43  and  44  show,  again  for  the  same  data  from  Figs.  41  and  42,  the  number  of  camera  views 
and  the  camera-based  precision  rmscam,  calculated  as  shown  in  Section  6.9.  (Recall  that  rmscam  is  a 
measure  of  the  disagreement  between  multiple  camera  views  for  a  given  point  on  the  test  article.) 

The  distribution  in  Fig.  44  is  essentially  the  same  as  that  in  Fig.  42,  though  spottier  (note  that  the 
rmscam  diagnostic  cannot  be  computed  for  points  that  are  viewed  by  only  a  single  camera).  The 
similarity  of  the  distributions  in  Figs.  42  and  44  indicates,  again,  that  the  main  source  of  error  in  the  data 
is  photon  shot  noise.  In  particular,  the  similarity  of  the  distributions  in  Figs.  42  and  44  indicates  that  there 
are  no  unexplained  view-angle  dependencies  in  the  illumination-excitation-emission-detection  chain  of 
processes  that  forms  the  basis  for  a  PSP  measurement. 
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Figure  41,  Mapped  Pressure  Distribution 
on  Space  Shuttle  Model 


Figure  42.  Pixel-Based  Precision  Associated 
with  Fig.  41 


Figure  43.  Number  of  Camera  Views  Figure  44.  Camera-Based  Precision 

for  Fig.  41  Associated  with  Fig,  41 

Overall,  the  ratios  between  the  three  diagnostie  measures  employed  was  found  to  be  on  the  order  of 
rmspix  :  rmscam  :  rmstap  =  3:4:8.  (57) 

This  result  makes  sense  qualitatively:  the  pixel-based  error  should  be  the  smallest  because  it  is 
essentially  a  measure  of  photon  shot  noise  and  is  least  affected  by  systematic  errors  such  as  view  factors, 
uncompensated  temperature  errors,  and  other  effects.  On  the  other  hand,  the  comparison  between  paint 
values  and  tap  values  captures  all  random  and  systematic  errors  in  the  measurement  scheme. 

As  to  why  the  rmstap  measure  is  three  times  as  large  as  the  rmspix  value,  the  reason  for  this  is  not 
clear  at  this  point,  but  it  is  suspected  that  this  ratio  is  uncommonly  high  because  of  the  geometrical 
complexity  of  the  test  article  and  the  fact  that  many  taps  are  located  in  or  near  comers,  where  it  is  not  easy 
to  assign  a  PSP-based  pressure  value.  Indeed,  if  the  rmstap  diagnostic  is  calculated  for  the  best  75 
percent  of  all  data  (i.e.,  the  PSP  data  that  best  match  the  available  taps),  the  new  value  of  rmstap  is  on 
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the  order  of  rmscam.  This  is  a  more  signifieant  improvement  in  rmspix  than  would  be  expected  if  the 
PSP-tap  differences  were  normally  distributed  (which  would  leave  a  new  ratio  of  5  :  8  vs.  4  :  8  for 
rmstap  :  rmscam).  In  other  words,  the  PSP-tap  differences  are  not  normally  distributed. 

Absolute  values  for  rmstap  are  quoted  in  Ref  66  (see,  e.g..  Fig.  9  in  Ref  66).  These  values  range 
from  0.5  kPa  to  1.3  kPa,  depending  on  test  conditions. 

Finally,  Fig.  45  shows  a  comparison  between  the  old  method  of  processing  data  (as  described  in 
Section  6.3)  and  the  new  method  from  Section  6.4.  This  comparison  shows  that  a  seam  along  the 
boundary  line  between  adjacent  camera  views  in  the  old  scheme  is  no  longer  visible  in  the  new  result. 
This  is  because  the  new  processing  scheme  employs  a  weighted  camera  average  rather  than  selection  of  a 
best  view,  and  because  per-camera  calibration  offsets  are  no  longer  required  since  converting  signal  ratios 
to  pressures  does  not  take  place  until  after  the  ratios  have  been  mapped  onto  the  3D  grid. 


Figure  45.  Example  of  Elimination  of  Seam  Resulting  from  Improved  Processing 
6.6  DISCUSSION 

It  has  been  shown  that  AEDC’s  new  lifetime-based  PSP  system  is  capable  of  producing  high-quality 
test  data  and  that  these  data  can  be  delivered  to  the  test  customer  in  near  real  time  during  the  PSP  test.  It 
has  also  been  shown  that  a  number  of  improvements  in  processing  methodology  have  still  been  possible, 
with  regard  to  both  quality  and  speed.  These  changes  involve  a  modification  of  the  order  in  which  the 
principal  processing  steps  are  executed,  improved  mapping  of  2D  image  data  to  the  3D  grid,  and 
averaging  data  from  multiple  camera  views  as  opposed  to  selecting  a  best  view.  In  particular,  converting 
signal  ratios  to  pressure  only  after  the  data  have  been  mapped  to  the  3D  grid  has  been  shown  to  both 
simplify  the  analysis  process  and  improve  the  final  results  (for  example,  the  elimination  of  seams  in  the 
processed  data  where  adjacent  camera  views  overlap). 

It  has  also  been  shown  that  additional  diagnostic  information  (beyond  the  conventional  comparison  of 
PSP  results  to  pressure  taps)  may  be  extracted  as  part  of  the  image-processing  task.  In  particular,  two 
diagnostic  measures,  rmspix  and  rmscam  (which  measure  the  precision  of  the  processed  data  on  a  per- 
pixel  and  a  per-camera  basis,  respectively),  have  been  implemented.  These  statistics  confirm  that  the 
precision  of  the  processed  data  for  the  current  lifetime-based  PSP  system  is  governed  essentially  by 
photon  shot  noise.  (This  statement  would  not  be  true  for  the  old  intensity-based  system.)  The  important 
implication  of  this  observation  is  that  it  should  be  possible  to  further  improve  the  precision  and  accuracy 
of  the  lifetime-based  PSP  system  by  increasing  the  radiant  flux  of  the  excitation  light.  Indeed,  plans  are  in 
progress  to  replace  the  current  LED  banks  with  higher  brightness  LED’s  that  have  recently  become 
commercially  available. 
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6.8  APPENDIX  6A:  IMPROVED  MAPPING  SCHEME 

In  Fig.  39,  assume  that  the  pixel  values  in  the  image  [representing  the  signal  ratio  from  Eq.  (55)]  are 
given  by  Ry,  where  i  and  j  are  indices  along  the  axes  of  the  image.  Furthermore,  assume  that  the  (floating 
point)  coordinates  of  the  projected  grid  point  g  are  given  by  {Ug,  Vg)  and  that  the  coordinates  of  the  centers 
of  the  pixels  in  the  vicinity  of  {ug,  Vg)  are  given  by  (m„  Vj).  Assignment  of  a  mapped  value  to  grid  point  g 
may  now  be  thought  of  as  an  interpolation  scheme  among  the  surrounding  pixel  values. 

To  first  approximation,  the  image  values  R(u,  v)  for  a  set  of  coordinates  (u,  v)  in  the  vicinity  of  (ug,  Vg) 
may  be  approximated  by  a  bilinear  slope  of  the  form 

R(u,v)  =  a  +  b(u -Ug)  +  c(v -v^)  ,  (58) 


where  a,  b,  and  c  are  coefficients  to  be  obtained  by  performing  a  least-squares  fit  of  Eq.  (58)  to  the 
measured  values,  Ry,  at  the  locations  (m„  Vy).  Equation  (58)  is  structured  in  such  a  way  that  the  coefficient 
a  yields  the  desired  mapped  value  with,  in  principle,  subpixel  resolution. 

Solution  of  the  least-squares  problem  leads  to  the  linear  system 


^  a] 

b 

= 

l^c  J 

(59) 


where  S\  El,  E(vy-Vg),  E(i//-i/g)  ,  E(i//-i/g)(vy-Vg),  Syy  E(vy-Vg)  ,  Sf  ERy, 

=  l.(ui-Ug)Rij,  and  Syy  =  E(vy-Vg)Ry.  All  of  the  sums  extend  over  the  same  set  of  pixels  (for  example,  the 
shaded  set  of  pixels  in  Fig.  39)  to  obtain  the  mapped  value  (i.e.,  the  value  a)  at  grid  point  g. 


The  unbiased  sample  variance  of  the  bilinear  fit  may  be  calculated  as 


N-3 


X  [;;(»„  V, )  -  A,  ]  =  =  ^  (S„  -  aS,  -  bS,„  -  cS^ ) , 


(60) 


where  Syr  =  '^Ri/  and  N  =  Si  =  HI  (i.e.,  the  number  of  pixels  in  the  sum).  The  factor  N-3  reflects  the  fact 
that  there  are  three  fit  coefficients  in  Eq.  (58). 

In  conjunction  with  Eq.  (56),  the  variance  from  Eq.  (60)  allows  the  rmspix  diagnostic  from  Fig.  38 
(for  grid  point  g)  to  be  calculated  as 


(rmspix)g 


k=l 


a 


g  * 


(61) 


For  best  results,  the  number  of  pixels  A  should  not  be  too  small.  (For  the  2048  by  2048  images  used 
in  this  work,  5  by  5  regions  of  pixels  were  found  to  give  good  results.)  For  completeness  it  should  be 
noted  that  other  types  of  mapping  are  possible,  for  example,  texture  mapping.^^ 


6.9  APPENDIX  6B:  OVERLAPPING  CAMERA  VIEWS 

In  the  old  mapping  scheme  from  Section  6.3,  the  best  view  for  a  grid  point  g  from  among  multiple 
camera  views  c  was  selected  on  the  basis  of  the  cosine  of  the  viewing  angle,  given  by 
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^cg=(«c*«g)-  (62) 

Here  ric  is  a  normal  vector  along  the  line  of  sight  for  camera  c,  and  rig  is  the  surface  normal  at  grid 
point  g.  In  the  old  scheme,  the  view  with  the  largest  value  of  Wcg  was  selected.  In  the  new  scheme  from 
Section  6.4,  a  weighted  average  is  calculated  in  the  case  of  multiple  camera  views,  according  to 

(63) 

C  I  C 

where  the  sum  is  over  all  camera  views  c  that  can  include  grid  point  g,  and  Reg  is  the  signal  ratio  for  grid 
point  g  as  measured  by  camera  c. 

The  degree  to  which  the  individual  ratios  Reg  vary  from  the  calculated  average  is  another  measure  of 
the  variance,  already  expressed  in  Eq.(60).  In  the  present  case,  this  variance  may  be  computed  according 
to 

y  z  ’  (64) 

where  Ne  is  the  number  of  camera  views.  The  variance  from  Eq.(64)  may  be  substituted  into  Eq.(61),  this 
time  to  obtain  the  camera-based  diagnostic  rmscam  from  Fig.  38. 
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7,0  ASSIMILATION  OF  PHYSICAL  CHEMISTRY  MODELS  FOR  LIFETIME 
ANALYSIS  OF  PRESSURE-SENSITIVE  PAINT^^ 


7.1  OVERVIEW 

Lifetime  analysis  teehniques  from  the  physieal  ehemistry  literature  are  assimilated  for  use  with 
lifetime-based  PSP  measurements.  Effects  of  finite  pulse  duration  and  residual  fluorescence,  typically 
avoided  in  physical  chemistry,  are  included  through  the  use  of  an  Excitation  Response  Function.  Three 
models  are  reviewed  and  applied  to  the  fluorescence  decay  of  PtTFPP  in  FIB,  which  is  the  current 
industry  standard  for  large-scale  PSP  testing  in  the  U.S.  These  are  the  Discrete  Exponential  Model,  the 
Forster  Decay  Model,  and  the  Maximum  Entropy  Model.  It  is  shown  how  these  models  invite  the 
reconstruction  of  the  decay  rate  distribution  of  the  paint  from  time-resolved  measurements  of  the 
fluorescence  response.  Results  of  the  analysis  suggest  that  the  decay-rate  distribution  (and  hence  the 
lifetime  spectrum)  of  the  paint  is  bimodal,  with  each  component  obeying  a  linear  Stem-Volmer 
dependence. 

7.2  INTRODUCTION 

Over  the  last  decade,  significant  interest  has  developed  in  the  use  of  luminescent  paints  for 
aerodynamic  testing. This  is  particularly  true  for  the  use  of  pressure-sensitive  paints  (PSPs),  which  rely 
on  oxygen-quenching  of  fluorescence  to  produce  a  pressure-dependent  luminescence  signal.  To  lowest 
order,  this  process  is  described  by  a  Stem-Volmer  relation  that  expresses  the  fluorescence  decay  rate,  y,  in 
terms  of  the  pressure,  P,  as 

y  =  -  =  — -tKsv^.  (65) 

X  To 

Here  x  is  the  lifetime  in  air,  Xo  is  the  lifetime  in  vacuum,  and  Ksv  is  the  bimolecular  quenching 
constant  in  air.  Higher  order  terms  are  frequently  included  in  Eq.  (65),  though  their  origin  is  usually  not 
well  understood.™’’^ 

In  most  PSP  measurements  neither  the  decay  rate,  y,  nor  the  lifetime,  x,  is  measured  directly.  Instead, 
a  ratio  is  formed  between  images  at  “wind-off’  and  “wind-on”  conditions  when  continuous  illumination  is 
used,  or  between  images  acquired  at  separate  gates  following  the  illumination  pulse  in  the  case  of  pulsed 
illumination.  These  two  techniques  are  referred  to  as  intensity-based  PSP  and  lifetime-based  PSP, 
respectively.  In  both  cases,  the  signal  ratio  is  expressed  as  a  function  of  pressure  by  means  of  a  Stern- 
Volmer-type  relation.  To  use  this  relationship  as  the  basis  for  a  pressure  measurement,  one  does  not 
necessarily  have  to  have  a  thorough  understanding  of  the  nature  of  the  relationship  (for  example,  linear 
vs.  nonlinear);  as  long  as  a  suitable  calibration  can  be  performed.  This  approach  has  been  fhiitful  since 
the  early  days  of  PSP  measurements.™'™ 

A  longstanding  complication  in  PSP  measurements  is  that  the  signal  ratios  depend  not  only  on 
pressure,  but  also,  to  some  extent,  on  temperature.*’’^'™  Recent  work  on  lifetime-based  PSP  suggests  that 
this  unwanted  temperature  dependence  may  be  turned  to  an  advantage  by  extending  the  standard  two-gate 
measurement  scheme  to  one  that  employs  three  or  four  gates.  This  allows  two  signal  ratios  to  be 

formed.  Under  suitable  conditions,  this  allows  measurement  of  both  pressure  and  temperature 
(see  Section  9.0). 


^^Previously  published  as  “Assimilation  of  Physical  Chemistry  Models  for  Lifetime  Analysis  of  Pressure-Sensitive 
Paint,”  AlAA  Journal,  Vol.  43,  No.  3,  March  2005,  pp.  465-471.  An  earlier  (and  longer)  version  of  this  paper  was 
presented  as  AIAA  Paper  2004-0880,  42nd  AIAA  Aerospace  Sciences  Meeting  and  Exhibit,  Reno,  NV,  January  5-8, 
2004. 
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Questions  remain  on  how  to  bring  such  a  P-T  measurement  scheme  to  fruition.  However,  it  is  clear 
that  the  fluorescence  response  of  the  paint  must  not  be  single-exponential:  multiple  measurements  on  a 
single-exponential  decay  would  merely  recover  the  same  y  or  t  in  different  ways,  leaving  the  P-T 
ambiguity  unresolved.  As  it  turns  out,  nonsingle-exponential  fluorescence  decay  of  PSPs  appears  to  be 
the  rule,  rather  than  the  exception  (see  Refs.  12,  13,  75,  and  78  and  Section  7.4  of  this  report).  Even  so, 
most  PSP  analysis  to  date  has  been  based  on  single-exponential  behavior.^’^*’’^  This  approach  falls  short  if 
the  optimization  of  a  multigate,  lifetime -based  measurement  scheme  is  to  be  approached  from 
fundamental  principles.  Instead,  it  is  desirable  to  start  out  from  the  known  lifetime  characteristics  of  the 
paint  or  (in  the  case  of  paint  development)  from  the  desired  paint  characteristics.  In  particular,  it  is 
desirable  that  these  lifetime  characteristics  be  cast  in  the  form  of  a  mathematical  model  that  can  serve  as 
the  basis  for  optimization  of  the  operational  parameters  of  a  lifetime -based  PSP  system,  for  derivation  of 
suitable  calibration  equations,  and  for  analysis  of  PSP  measurements  in  general.  Such  lifetime  analysis 
models  have  been  in  use  for  some  time  in  physical  chemistry  and  associated  fields  (chemical  physics, 
analytical  chemistry,  sol-gel  physics,  and  biophysics).  The  purpose  of  Section  7.0  is  to  assimilate  some  of 
these  models  for  use  with  PSP  (Section  7.3)  and  (in  Section  7.4)  to  demonstrate  the  use  of  these  models 
for  PtTFPP  in  FIB  [Pt(II)  meso-tetra(Pentafluorophenyl)porphine  in  fluoro-isopropyl-butyl],  which  is  the 
current  industry  standard  for  PSP  testing  in  the 


7,3  LIFETIME  MODELS  FROM  PHYSICAL  CHEMISTRY 

In  Section  7.3.1,  the  notion  of  paint  as  a  microheterogeneous  system  is  introduced,  leading  to  a 
general  formulation  for  the  response  of  such  a  system  to  delta-fiinction  excitation.  In  Sections  7.3.2,  7.3.3, 
and  7.3.4,  three  specific  realizations  of  this  general  model  are  described.  Section  7.3.5  deals  with  the 
situation  in  which  the  width  and  repetition  rate  of  the  excitation  pulse  can  not  be  neglected.  Section  7.3.6 
discusses  numerical  aspects  of  the  reconstruction  of  the  lifetime  spectrum  of  the  paint  from  a  time- 
resolved  measurement  of  the  fluorescence  signal. 


7,3,1  Paint  as  a  Microheterogeneous  System 

A  central  notion  in  the  description  of  quenched  fluorophore  systems  is  that  not  all  fluorophore 
molecules  experience  the  same  interaction  with  the  quenching  species  (in  the  case  of  current  PSPs,  this 
would  be  oxygen)  and  that  this  heterogeneity  exhibits  itself  on  a  microscopic  scale. Such  a  system 
may  be  described  by  partitioning  the  fluorophore  density,  Wf,  with  respect  to  the  resulting  distribution  of 
decay  rates,  y,  as 


(66) 


That  is,  Uf  is  replaced  by  a  distribution  fimction  n^'{y)  that  has  Uf  as  its  area.  Equivalently,  this 
partitioning  may  be  described  with  respect  to  the  lifetime,  x,  where  x  =  l/y.*°  The  fluorescence  response  of 
such  a  system  to  delta-function  excitation  may  be  expressed  as  the  Faplace  transform 


55(0  =  jQ”«y(Y)e 


(67) 


where  «y(y)  is  the  distribution  function  ny'(y)  from  Eq.  (66),  multiplied  by  a  constant: 

ny{T)  =  Cn\{y)  .  (68) 

The  constant  C  may  be  written  (see  Section  7.6)  as  the  product  of  the  area  viewed  by  the  detector,  the 
effective  thickness  of  the  sample,  the  excitation  photon  flux  at  the  surface,  the  optical  cross  section  for 
excitation,  the  quantum  efficiency  for  population  of  the  fluorescing  state,  the  rate  of  fluorescence  decay, 
the  efficiency  of  detection  of  emitted  fluorescence  photons,  the  detection  time,  and  the  number  of 
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excitation  pulses  over  which  the  signal  is  integrated.*®  Implicit  in  Eq.  (68)  is  the  assumption  that  the 
excitation  process  is  linear  (see  Section  7.7). 

From  Eqs.  (66)  through  (68),  it  follows  that  the  area  under  the  distribution  nfy)  is  given  by 

J  =  j”n,(Y)dy  =  55(0)  =  Cnf.  (69) 

This  area  should  be  independent  of  quenching,  a  principle  that  might  be  paraphrased  as  “conservation 
of  implied  fluorophore  density.” 

7.3.2  Discrete  Exponential  Model 

When  a  single  exponential  (with  to-be-determined  amplitude  and  decay  rate)  fails  to  produce  a 
satisfactory  fit  to  fluorescence  decay  data,  it  is  common  practice  to  add  a  second,  third,  and  (in  rare  cases) 
a  fourth  term.'^’^'’^^  This  model  is  known  as  the  Discrete  Exponential  Model  (DEM)  and  may  be 
expressed  as 

K 

5'5,  DEM  (0  =  Z  exp(-YA;0 .  (70) 

k=l 

The  number  of  components,  K,  may  be  assigned  on  the  basis  of  an  a  priori  model  (for  example,  when 
the  components  are  associated  with  different  fluorophore  molecules),  or  chosen  empirically,  to  get  a 
satisfactory  fit  to  the  data.  Figure  46  illustrates  the  associated  decay  rate  distribution  for  the  case  K  =  3. 
The  “area”  under  this  distribution  is  given  by  the  relation 

K 

-^DEM^Z'^A;-  C71) 

k=\ 


«Y(y) 


Figure  46.  Illustration  of  DEM  Decay  Rate  Distribution 
7.3.3  Forster  Decay  Model 

If  the  range  of  interaction  distances  between  fluorophore  molecules  and  quenchers  is  given  by  a 
power-law  distribution,  the  delta-function  response  of  the  system  may  be  written  as*^'®® 

-^e.FDMCO  = '^'o  exp(-at-^t^^^),  (72) 

where  Sq  is  the  peak  signal,  a  is  the  Stem-Volmer  quenching  rate  from  Eq.  (65)  with  a  =  y  =  1/x,  and  ^  is  a 
measure  of  the  heterogeneity  of  the  quenching  process.  The  exponent  1/2  in  Eq.  (72)  results  for  a  three- 
dimensional  environment.  Exponents  1/3  and  1/6  apply  to  two  dimensions  and  one  dimension, 
respectively.®®  Sometimes  the  exponent  itself  is  used  as  a  fitting  coefficient.®'  Of  particular  interest  to  PSP 
research  is  work  by  Draxler  et  al.,**’*®  who  found  that  the  oxygen  quenching  of  Ru(dpp)  in  a  series  of 
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polymers  (i.e.,  a  PSP-like  system)  could  be  described  better  by  the  Forster  Decay  Model  (FDM)  from  Eq. 
(72)  than  by  the  DEM  from  Eq.  (70). 

Taking  the  inverse  Eaplace  transform  of  Eq.  (72),  it  follows  (see  Ref  92)  that  the  associated  Forster 
decay  rate  distribution  is  given  by: 


«y(Y) 


4^0  J  b 

4nb^  [lyjj-a  ^ 


(73) 


where  G{x)  =  exp(-x^),  and  ny(y)  =  0  for  y  <  a.  Figure  47  illustrates  this  distribution,  which  peaks  at  y  = 

a  +  b^l6  and  has  an  approximate  height  OSlSS^Ib^  and  an  approximate  full  width  at  half  maximum  given 
by  0.45^^.  The  area  under  this  distribution,  defined  by  Eq.  (69),  is  precisely  equal  to  the  amplitude  Sq,  i.e., 

^FDM  =  Sq. 


FDM 


7.3.4  Maximum  Entropy  Model 

A  model  that  is  formally  equivalent  to  the  DEM  is  one  in  which  the  decay  rate  distribution,  ny(y),  is 
represented  by  a  histogram  of  preselected  decay  rates,  y^,  with  bin  widths  Ay^  and  to-be-determined  bin 
areas  (see  Fig.  48). 


Figure  48.  Illustration  of  MEM  Decay  Rate  Distribution 

The  delta-function  response  of  such  a  system  can  be  written  as  [cf  Eq.  (70)]: 

M 

5'8,MEm(0  =  expC-TnO.  (74) 

H=1 
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where  M  is  the  number  of  bins  in  the  histogram.  Typical  values  of  M  are  in  the  range  from  100  to  200. 
Fitting  Eq.  (74)  to  a  measured  fluorescence  response  is  a  linear  problem  (because  the  values  are 
preselected),  but  one  that  is  highly  ill  conditioned.*"^’^^'^"'  To  obtain  bin  areas,  that  are  physically 
realistic  (not  negative  and  not  strongly  oscillatory),  it  is  necessary  to  impose  a  constraint.  In  the 
Maximum  Entropy  Model  (MEM),  this  constraint  is  taken  to  be  the  Shannon-Jaynes  entropy^^''®° 

M 

E  =  A) ln(a^  /  A)  ,  (75) 

^l=l 

where  A  is  the  actual  or  estimated  area  under  the  distribution,  given  by  [cf  Eq.  (71)  for  the  DEM]: 

M 

-^MEM  = 

H=1 

The  discretized  values  of  the  distribution  function  ^^(y)  are  given,  in  terms  of  the  bin  areas  and  the 
bin  widths  Ay^,  by  =  a^/Ay^. 

7.3.5  Excitation  Response  Function 

In  physical  chemistry,  lifetime  studies  of  quenched  fluorophore  systems  typically  are  performed  using 
short-pulsed  lasers  so  that  the  assumption  of  delta-function  excitation  is  satisfied.  However,  in  PSP 
studies  this  may  not  be  the  case  (for  example,  when  pulsed  EEDs  are  used  as  the  excitation  source). 
Moreover,  to  minimize  PSP  measurement  time,  the  optimum  pulse  repetition  frequency  may  be  such  that 
small  amounts  of  residual  fluorescence  (from  prior  pulses)  are  tolerated.  Both  of  these  complications  can 
be  accounted  for  by  convolving  the  presumed  delta- function  response  of  the  paint  [from  Eq.  (67),  (70), 
(72),  or  (74)]  by  the  normalized  excitation  pulse  shape,  p{t),  to  arrive  at  the  measured  signal,  S{t), 
according  to 


S(0  =  j_'^M^o)^8(^-^o)d?o-  (77) 

The  lower  limit  of  the  integral  is  set  to  -oo  to  indicate  that  the  measured  fluorescence  response  at  time 
t  extends  over  the  entire  excitation  history  of  the  system,  including  residual  fluorescence  from  prior 
pulses. 

By  combining  Eqs.  (67)  and  (77),  one  can  write  the  signal  S{t)  as 

C  00 

‘^(0  =  |o  ^(Y)'^y(^;Y)dY,  (78) 

where  5y(t;y)  is  the  Excitation  Response  Function  (ERF)  for  monoexponential  decay  at  a  rate  y,  given  by 

Sy{f,y)  =  .  (79) 

The  ERF  may  be  calculated  analytically  (see  Section  7.8)  for  simple  pulse  shapes,  pulse  trains,  and 
waveforms,  including  the  case  of  continuous  excitation,  with  p{t)  =  1  and  Sy(t;y)  =  l/y.*°  Alternatively,  the 
ERF  must  be  calculated  by  numerical  convolution  if  the  pulse  shape  p{i)  is  measured  experimentally  and 
cannot  be  described  adequately  by  a  simple  analytical  expression.^' 

7,3,6  Reconstruction  of  Decay  Rate  Distribution  from  Fluorescence  Response 

From  the  foregoing  it  follows  that  fitting  a  measured  fluorescence  response  can  be  interpreted  as 
reconstructing  the  underlying  distribution  of  fluorescence  decay  rates.  To  do  so,  assume  that  the 
fluorescence  response  S{t)  has  been  sampled  to  yield  a  set  of  N  values  for  A,  and  that  the  pulse  shape  p{t) 
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has  been  measured,  with  the  result  that  the  ERF  from  Eq.  (79)  can  be  calculated  for  any  desired  value  of 
y.  Fitting  the  above  lifetime  models  to  the  fluorescence  response  data  can  then  be  cast  as  a  least-squares 
fitting  problem,  in  which  the  quantity 

(80) 

/=1 

is  to  be  minimized.  Here  represents  the  fitted  signal  value  at  point  /,  and  a,  is  the  estimated 
uncertainty  for  this  point,  proportional  to  the  square  root  of  Si  if  photon  shot  noise  is  assumed. 

In  the  case  of  the  DEM  and  EDM,  standard  nonlinear  regression  can  be  used  to  determine  optimum 
values  of  the  fit  parameters.'®'  In  the  DEM  from  Eq.  (70),  these  are  the  K  amplitudes,  ak,  and  associated 
decay  rates,  y^;  in  the  EDM,  these  are  the  parameters  5*0,  a,  and  b  from  Eq.  (72).  In  the  case  of  the  MEM, 
the  to-be-minimized  quantity  is  where  is  given  by  Eq.  (80),  E  is  given  by  Eq.  (75),  and  is  a 

to-be-determined  Eagrange  multiplier.  Because  the  number  of  to-be-determined  bin  areas,  Uk,  in  Eq.  (74) 
typically  is  in  the  range  from  100  to  200,  a  special  algorithm  is  required.®^’'®®  [The  original  Maximum 
Entropy  Method  (MEM)  algorithm  was  developed  for  image  reconstruction  in  astrophysics'®^  and  has 
since  been  adapted  for  lifetime  analysis.®^] 

Years  of  research  in  physical  chemistry  (beginning  with  a  seminal  paper  by  James  and  Ware'®^)  have 
made  it  clear  that  it  can  be  difficult  to  select,  on  the  basis  of  the  magnitude  of  fitted  residuals  alone,  the 
best  model  for  a  particular  system.  For  example,  even  if  the  true  dynamics  are  given  by  a  Forster  model,  a 
two-component  DEM  can  give  a  good  fit  to  the  data.*®'**  For  this  reason,  it  is  desirable  to  consider  not 
only  the  magnitude  of  the  residuals,  but  also  the  behavior  of  the  fitted  parameters  as  a  function  of 
quencher  density.  This  is  illustrated  in  Section  7.4. 

7,4  LIFETIME  ANALYSIS  OF  PtTFPP  IN  FIB 

In  this  section,  the  theory  from  Section  7.3  is  applied  to  fluorescence  decay  measurements  on  the 
pressure-sensitive  paint  PtTFPP  in  FIB,  which  was  developed  at  the  University  of  Washington®  and  is 
marketed  by  Innovative  Scientific  Solutions,  Inc. 

7,4,1  Experiment 

Painted  aluminum  coupons  were  mounted  in  a  calibration  chamber  in  which  air  pressure  and  sample 
temperature  were  precisely  controlled.  Shop  air  with  a  specified  dewpoint  of  -50°F  was  used  to  avoid 
complications  caused  by  humidity  effects.  Samples  were  illuminated  at  a  45-deg  angle  by  a  commercial 
array  of  blue  EEDs  with  peak  emission  at  465  nm,  resulting  in  excitation  of  the  PtTFPP  Soret  band.  The 
illumination  level  was  held  constant  to  better  than  one  percent  over  the  course  of  the  measurements. 
Fluorescence  from  the  triplet  state,  filtered  at  650  nm  using  a  40-nm  FWHM  custom  bandpass  filter  with 
70-percent  transmission,  was  measured  by  a  fast  photomultiplier  tube  (PMT).  The  observation  direction 
was  also  at  a  45-deg  angle  to  the  sample  normals,  opposite  from  the  incident  light  direction. 

The  FED  source  was  pulsed  at  4.95  kHz,  giving  a  train  of  nearly  square  light  pulses  spaced  at  202  ps, 
each  with  a  duration  of  approximately  10  ps.  Output  from  the  PMT  was  sampled  at  0.1 -ps  intervals  from 
-0.5  ps  to  about  -i-lOO  ps  relative  to  the  trigger  input,  resulting  in  1002  samples  per  curve.  Curves  were 
summed  over  1 00  pulses  using  a  digital  oscilloscope.  Preprocessing  of  the  data  was  performed  to  subtract 
a  baseline  value  and  to  determine  the  precise  start  of  the  optical  pulse.  Calibrations  were  carried  out  over 
a  matrix  of  25  pressures  (from  20  to  2100  psf),  each  at  10  temperatures  in  the  range  from  40  to  120°F. 

Figure  49  shows  the  measured  optical  pulse  shape,  p{t),  which  has  a  rise  time  of  about  0.5  ps  and  is 
followed  by  a  small-amplitude  tail  (0.6  percent  of  peak  signal)  that  extends  about  13  ps  beyond  the  main 
cutoff 
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Figure  50  shows  the  measured  fluorescence  signal,  S(t),  for  several  pressures  at  room  temperature. 
Data  are  plotted  on  both  a  linear  and  a  log  scale.  From  the  decay  part  of  the  fluorescence  in  Fig.  50b,  it  is 
evident  that  the  fluorescence  response  is  not  single-exponential.  Close  examination  of  the  data  at  the 
lowest  pressures  (i.e.,  the  data  with  the  longest  lifetimes)  also  reveals  the  effect  of  residual  fluorescence 
from  pulses  prior  to  the  trigger  pulse,  with  a  maximum  relative  contribution  of  four  percent  of  the  peak 
signal  (see  Fig.  3  in  Ref.  81). 
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b.  Log  Scale  Decay  Portion  Only 
Figure  50.  Selected  Fluorescence  Response  Curves  for  PtTFPP  in  FIB 

7.4.2  Data  Analysis 

All  analysis  was  based  on  least-squares  fitting  of  the  measured  fluorescence  response  curves, 
including  the  excitation  phase  (See  Fig.  51).  Numerical  convolution  was  used  to  take  into  account  the 
measured  pulse  shape  and  pulse  repetition  rate.*'  Shot  noise  was  assumed,  with  relative  weights  a, 
proportional  to  the  square  root  of  the  signals  5,.  A  Levenberg-Marquardt  algorithm  (with  explicit 
calculation  of  derivatives  with  respect  to  the  fit  parameters)  was  used  for  the  DEM  and  FDM  models.  The 
MEM  algorithm  was  adapted  from  Ref  100,  with  some  insights  gleaned  from  Refs.  93  and  98. 
Convergence  times  for  all  three  models  were  typically  on  the  order  of  several  seconds  per  data  set  on  a 
750-MFIz  Pentium-based  processor.  A  simple  peak-finding  algorithm  was  used  to  identify  peaks  in  the 
MEM  distributions  and  to  calculate  the  partial  area  fractions  under  the  peaks  relative  to  the  total  area 
under  the  distribution. 
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Time,  |is 

Figure  51.  Fit  Residuals  for  Data  from  Fig.  50  for  Four  Models 

Of  the  three  models  tested,  MEM  eonsistently  produeed  the  best  fit  to  the  data.  This  is  to  be  expeeted, 
beeause  of  the  large  number  of  adjustable  parameters.  Sample  fit  residuals,  as  a  pereentage  of  the  peak 
signal  (i.e.,  the  signal  at  the  end  of  the  exeitation  pulse),  are  shown  in  Fig.  51.  Insets  also  state  the  average 
RMS  fit  error,  s,  as  a  pereentage  of  the  peak  signal.  Clearly,  monoexponential  deeay  (“DEM-1”)  fails  to 
fit  the  data.  Differenees  between  the  two-eomponent  DEM  (“DEM-2”),  EDM,  and  MEM  residuals  are 
more  subtle,  with  the  MEM  residuals  showing  the  least  strueture. 

Figure  52  shows  the  relative  RMS  error,  s,  versus  pressure.  As  stated,  MEM  eonsistently  produees  the 
smallest  fit  errors.  Interestingly,  for  most  pressures  the  three-parameter  EDM  fits  produee  lower  residuals 
than  the  four-parameter  DEM  fits. 

Representative  deeay  rate  distributions  for  the  three  models  are  shown  in  Fig.  53.  The  vertieal  axis  of 
Fig.  53  has  two  sets  of  units:  one  for  the  DEM  amplitudes  and  one  for  the  EDM  and  MEM  distributions. 
In  the  EDM  ease,  the  spectrum  was  obtained  from  Eq.  (73)  after  the  data  were  fitted  using  Eq.  (72).  While 
the  EDM  spectrum  has  only  one  peak,  the  MEM  and  (of  course)  two-component  DEM  spectra  have  two. 
For  some  data  points,  a  small  third  peak  was  found  in  the  MEM  spectrum  at  small  values  of  y.^’  More 
accurate  measurements  will  be  required  to  confirm  that  this  third  peak  is  not  an  artifact  of  the 
measurement  and  analysis  schemes.  Eikewise,  a  three-component  DEM  fit  produced  only  a  marginal 
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Figure  52,  Rms  Fit  Residuals  as  a  Function  of  Pressure 


■y(=1/t],  n.s-1 

Figure  53.  Implied  Decay  Rate  Distributions  for  Data  from  Fig.  50 


improvement  over  the  two-component  DEM,  with  no  consistent  placement  of  the  third  peak  as  a  function 
of  pressure.  The  two-peaked  nature  of  the  MEM  spectrum,  however,  was  consistent  across  all  data  points, 
even  in  cases  in  which  the  two  peaks  were  not  fully  resolved.^’ 

Figure  54  shows  the  areas.  A,  under  the  ny(y)  distributions  as  a  function  of  pressure  for  the  three 
models.  As  argued  following  Eq.  (69),  these  areas  should,  ideally,  be  independent  of  pressure,  a  reflection 
of  the  fact  that  the  number  of  fluorophore  molecules  that  participate  in  the  excitation  process  is  constant, 
regardless  of  quenching.  Mild  variations  with  pressure  are  found  for  each  of  the  three  models. 


Figure  54,  Areas  Under  Distribution  Function  as  a  Function  of  Pressure 
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Figure  55  shows  the  FDM  parameters  a  and  b  from  Eq.  (72)  as  funetions  of  pressure.  Note  the 
tendency  of  a  towards  negative  values  at  large  pressures.  Figure  56  gives  both  the  MEM  and  DEM  peak 
positions  as  functions  of  pressure.  The  two  sets  of  peak  positions  are  in  reasonable  agreement,  both 
suggesting  a  linear  SV  dependence  as  expressed  by  Eq.  (65). 


Figure  55.  Forster  Parameters  a  and  6  as  a  Function  of  Pressure 


Figure  56.  Stern- Volmer  Diagram  with  Peak  Positions  vs.  Pressure 

Fitted  values  for  the  SV  intercepts  and  slopes  [expressed  in  the  form  of  the  vacuum  lifetimes,  to,  and 
bimolecular  quenching  constant,  Ksv,  from  Eq.  (65)]  are  listed  in  Table  9,  which  also  lists  the  estimated 
uncertainties  and  fitted  temperature  sensitivities  of  these  parameters.  Finally,  Fig.  57  gives  the  fractional 
contributions  of  the  two  peaks  in  the  DEM  and  MEM  distributions  to  the  total  area  under  the  distribution. 
As  is  the  case  in  Fig.  57,  there  is  good  agreement  between  the  DEM  and  MEM  results. 
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Table  9.  Fitted  Stern-Volmer  Parameters  for  PtTFPP  in  FIB  at  70°F 


Value 

Slope  (%/°F) 

Xo\  [IS 

51  ±4 

Unclear 

So 

21  ±5 

Unclear 

ps''Kpsr' 

0.035  ±0.004 

±0.11  ±0.04 

Kgy\  pS''KpSf^' 

0.168  ±0.015 

±  0.22  ±  0.02 

Figure  57.  Peak  Area  Fractions  Associated  with  Fig,  56 
7.4.3  Interpretation  of  Results 

Despite  the  fact  that  the  FDM  tends  to  produce  better  fits  to  the  data  than  does  the  DEM,  it  does  not 
appear  that  the  FDM  is  a  physically  plausible  model  for  the  fluorescence  decay  of  PtTFPP  in  FIB.  In 
particular,  the  dependence  of  the  SV  parameter  a  in  Fig.  55  is  opposite  the  rise  with  pressure  that  is 
expected,  and  was  found  by  Draxler  et  ah,  for  the  oxygen-quenching  of  Ru(dpp)  in  a  series  of  polymers.** 
Also,  unlike  what  was  found  by  Draxler  et  al.,  the  parameter  b  in  Fig.  55  is  strongly  dependent  on 
quencher  density  in  the  present  case. 

On  the  other  hand,  results  of  MEM  analysis  appear  to  confirm  that  a  two-component  DEM  captures 
the  essence  of  the  fluorescence  decay  process  correctly,  even  though  the  DEM  does  not  account  for 
nonzero  widths  of  the  peaks  in  the  decay  rate  distribution.  The  fact  that  the  MEM  produces  lower  fit 
residuals  than  does  the  DEM  indicates  that  the  widths  of  the  MEM  peaks  have  a  physical  basis,  even 
though  it  is  known  from  previous  work  on  MEM  that  noise  on  the  data  leads  to  broadening  of  the  peaks  in 
the  MEM  reconstruction.^®'^* 

What  remains  to  be  explained  is  why  the  fractional  areas  of  the  two  peaks  are  a  strong  function  of 
quencher  density.  This  behavior  has  previously  been  considered  nonphysical  from  the  perspective  of  a 
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microheterogeneous  system,  in  which  each  peak  is  associated  with  a  particular  group  of  fluorophore 
molecules.*^'**  An  alternative  explanation  is  suggested  in  Ref  104  (also  Section  8.0  of  the  present  report), 
which  involves  a  heretofore  unknown  intermediate  state  that  appears  to  play  a  role  in  the  fluorescence 
decay  process.  More  work  will  be  required  to  verify  the  existence  of  this  intermediate  state  (involving, 
perhaps,  detailed  spatial'®^  or  spectral  measurements’®^)  and  to  see  if  such  an  intermediate  state  plays  a 
similar  role  in  other  PSP  systems. 

7.5  DISCUSSION 

Three  lifetime  analysis  models  from  the  physical  chemistry  literature  have  been  reviewed  and  applied 
to  the  fluorescence  decay  of  PtTFPP  in  FIB,  which  is  the  current  industry  standard  for  large-scale  PSP 
testing  in  the  U.S.  Included  in  this  analysis  are  two  effects  that  typically  are  avoided  in  physical 
chemistry:  the  finite  width  of  the  excitation  pulse,  and  residual  fluorescence  from  prior  pulses.  Of  the 
three  models  studied,  the  Maximum  Entropy  Model  produces  the  best  fit  to  the  data  and  suggests  strongly 
that  there  are  two  components  in  the  fluorescence  decay.  The  associated  decay  rates  and  relative 
contribution  fractions  of  these  two  components  are  also  reproduced  by  a  Discrete  Exponential  Model, 
even  though  the  DEM  fit  residuals  are  larger  than  those  of  the  MEM.  A  Forster  Decay  Model  does  not 
appear  to  be  a  plausible  physical  model  for  PtTFPP  in  FIB,  even  though  it  produces  lower  fit  residuals 
than  does  the  two-component  DEM,  and  even  though  it  has  been  advocated  as  an  appropriate  model  for  a 
similar  system.** 

It  is  anticipated  that  insights  from  this  work  will  guide  the  further  development  of  lifetime-based  PSP 
techniques.  In  particular,  by  providing  tools  for  the  measurement  and  modeling  of  the  fluorescence  decay 
dynamics  of  PSP  systems,  it  should  become  possible  to  approach  the  optimization  of  lifetime-based  PSP 
measurements  in  a  fundamental  way.  This  should  be  of  particular  interest  to  current  efforts  in  which  the 
concept  of  a  two-gate  measurement  for  pressure  only  is  extended  to  one  that  encompasses  three  or  four 
gates,  while  exploiting  the  nonsingle-exponential  fluorescence  dynamics  of  the  paint  to  measure  both 
pressure  and  temperature. 

7.6  APPENDIX  7A:  DERIVATION  OF  PROPORTIONALITY  CONSTANT” 

The  constant  C  in  Eq.  (68)  may  be  derived  as  follows.  Consider  the  case  in  which  the  environment  of 
the  fluorophore  molecules  is  homogeneous,  in  the  sense  that  each  excited  fluorophore  molecule  has  the 
same  probability  of  deexcitation,  whether  it  be  through  emission  of  a  fluorescence  photon,  nonradiative 
decay,  or  quenching.  The  number  of  excited-state  fluorophore  molecules  in  the  observation  volume, 
N*(t),  may  then  be  expressed  by  the  rate  equation 

dN*(t)/dt  ^R(t)-rN*(t),  (81) 

where  R(t)  is  the  production  rate  of  excited-state  fluorophores  and  y  is  the  average  decay  rate  of  the 
excited  state.  The  solution  of  Eq.  (81),  for  the  initial  condition  A* (-oo)  =  0,  may  be  written  as 

=  (82) 

If  the  excitation  process  is  linear  (see  Section  7.7),  the  production  rate,  R(t),  may  be  written  in  terms 
of  the  normalized  pulse  shape,  p(t),  as 

R{t)  =  (83) 


^Material  in  Sections  7.6,  7.7,  and  7.8  is  taken  from  Ref.  23:  W.  Ruyten,  “Assimilation  of  Physical  Chemistry 
Models  for  Lifetime  Analysis  of  Pressure-Sensitive  Paint,”  AIAA  Paper  2004-0880,  42nd  AIAA  Aerospace 
Sciences  Meeting  and  Exhibit,  Reno,  NV,  January  5-8,  2004. 
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where  Ag  is  the  area  of  the  sample  viewed  by  the  detector,  ds  is  the  effective  thickness  of  the  sample  (see 
Section  7.7),  nf  is  the  ground-state  fluorophore  density,  Omax/^Vex  is  the  incident  photon  flux  at  the  surface 
(in  terms  of  the  power  density  Omax,  Planck’s  constant  h,  and  the  average  frequency,  Vex,  of  the  excitation 
light),  Gex  is  the  effective  excitation  cross  section,  and  Q  is  the  quantum  efficiency  for  conversion  of  the 
initially  excited  singlet  state  to  the  fluorescing  triplet  state.  The  ordering  of  the  factors  on  the  right-hand 
side  of  Eq.  (83)  follows  the  sequence  of  steps  that  a  fluorophore  molecule  experiences  during  its 
promotion  from  the  ground  state  to  the  fluorescing  state,  starting  with  tif,  the  number  of  available 
ground-state  molecules  in  the  observation  volume. 

Given  the  number  of  excited  fluorophores,  the  fluorescence  signal,  S(t),  may  be  written  as 

5'(0  =  A^*(0remhem^^^p^  (84) 

where  Tem  is  the  rate  of  fluorescence  (with  Tem  ^  l/x,  because  the  total  deexcitation  rate  involves 
nonradiative  decay  and  quenching  as  well),  pem  is  the  overall  efficiency  with  which  fluorescence  photons 
are  converted  to  units  of  signal  (e.g.,  mV),  At  is  the  sampling  interval  for  the  signal  S(t),  and  Np  is  the 
number  of  pulses  over  which  the  signal  is  integrated. 

Combining  Eqs.  (82)  through  (84)  and  integrating  over  the  decay  rate  distribution  ny(y)  from  Eq.  (67) 
then  yields  Eq.  (68),  with  a  constant  C  given  by 

C^A^d^  .  (85) 

^Vex 

7,7  APPENDIX  7B:  CONDITIONS  FOR  LINEAR  EXCITATION  AND  FLUORESCENCE 

The  following  assumptions  are  implicit  in  the  derivation  of  the  excitation  and  fluorescence  decay 
processes  from  Section  7.6. 

Assumption  1:  The  excitation  process  is  linear  with  respect  to  the  intensity  of  the  incident  light  (i.e.,  it 
is  not  a  multiphoton  process). 

Assumption  2:  The  decay  rate,  y,  is  independent  of  time.  In  the  context  of  PSP  measurements,  this  is 
to  say  that  the  measurement  takes  place  under  steady-state  conditions,  or  that  the  pressure  distribution 
varies  on  a  time  scale  that  is  much  longer  than  the  associated  lifetime  x  =  1/y. 

Assumption  3:  An  effective  sample  thickness,  ds,  can  be  assigned,  even  though  the  paint  may  be  (and 
typically  is)  optically  thick.  For  an  optically  thick  paint,  this  means  that  ds  is  less  than  the  actual  thickness 
of  the  paint. 

Assumption  4\  The  ground-state  fluorophore  density  that  is  available  for  excitation,  nf,  is  independent 
of  time.  This  assumption  is  satisfied  if,  and  only  if,  the  excitation  process  does  not  lead  to  a  significant 
depletion  of  the  ground-state  population.  In  practice,  this  sets  an  upper  limit  on  the  strength  of  the 
excitation  light,  under  which  a  linear  input-output  relation  holds  true. 

Assumption  5:  The  quantum  efficiency,  Q,  is  independent  of  time.  Specifically,  it  is  assumed  that 
conversion  from  the  initially  excited  singlet  state  to  the  fluorescing  triplet  state  is  instantaneous  compared 
to  the  lifetime  x. 

Assumption  6\  No  processes,  other  than  pumping  ground-state  atoms  via  the  singlet  state,  contribute 
to  the  production  of  excited-state  molecules.  Specifically,  this  assumes  that  delayed  fluorescence  is 
negligible,  by  which  triplet-state  molecules  undergo  an  internal  conversion  to  an  intermediate  state  and, 
after  some  delay,  are  converted  back  to  the  triplet  state. 

Assumption  7:  The  detection  efficiency,  pem,  is  independent  of  time.  Specifically,  this  assumes  that  a 
fast  detector  is  used,  whose  time  constant  is  fast  compared  to  the  lifetime  x. 
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Assumption  8:  Effective  values  of  Vex,  aex,  Q,  Eem,  and  rjem  may  be  assigned,  even  though  both  the 
excitation  process  and  the  emission  process  involve  wavelength  ranges  that  are  set  by  the  light  source,  the 
detector,  and  any  optical  fdters  used. 

All  these  assumptions  can  be  relaxed  in  favor  of  a  more  rigorous  treatment.  See,  for  example.  Ref. 
106  for  a  study  on  the  ability  of  a  PSP  measurement  to  resolve  a  dynamic  pressure  change,  or  Ref  107  for 
a  study  in  which  lifetime  analysis  of  a  PSP-like  system  is  performed  as  a  function  of  emission 
wavelength.  Also,  if  the  detector  response  is  not  instantaneous  but  given  in  terms  of  an  instrument 
response  function  (IRF),  this  IRE  can  be  convolved  with  the  ERF  to  obtain  a  total  response  function  for  a 
single  decay  rate  y. 

7.8  APPENDIX  7C:  SPECIAL  PULSE  SHAPES 

The  excitation  response  function  may  be  calculated  analytically  for  simple  pulse  shapes.  Table  10 
lists  the  results  for  excitation  by  a  single  delta  function  (with  effective  duration  Atp),  a  train  of  such  delta 
functions  spaced  at  a  cycle  time  4,  a  single  square  pulse  with  duration  tp,  a  rounded  square  pulse  with  the 
same  duration  but  with  a  rise-time  Xr,  a  sinusoid  with  frequency  co  =  2n/tc  and  modulation  depth  Mq  (or 
2mo),  and  a  continuous  wave. 


Table  10.  Excitation  Response  Functions  for  Selected  Pulse  Shapes 


Pulse  Shape,  p(t) 
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Figure  58  gives  a  graphical  illustration  of  the  pulse  shapes  p{t)  and  resulting  ERFs  s-,{t)  from  Table  10 
(i.e.,  the  ERFs  are  expressed  in  terms  of  the  lifetime  x  rather  than  the  decay  rate  y),  for  both  fast  and  slow 
decay.  The  relative  vertical  scaling  between p{t)  and  s-c(t)  in  Fig.  58  (i.e.,  between  quantities  with  different 
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units)  varies  by  row,  but  not  by  column.  Thus,  for  example,  in  Fig.  58c,  a  longer  lifetime  produces  not 
only  a  slower  rise  and  a  slower  decay,  but  also  a  higher  peak  signal.  Note  the  effects  of  residual 
fluorescence  in  Fig.  58b  for  the  long-lifetime  case.  Figure  58e  depicts  the  scenario  encountered  in 
frequency-domain  spectroscopy,  in  which  measurement  of  the  phase  shift  cpa,  as  a  function  of  the 
modulation  frequency  ®  forms  the  basis  for  lifetime  analysis. 


Figure  58.  Pulse  Shapes  and  ERFs  from  Table  10 
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8.0  OXYGEN  QUENCHING  OF  PTTFPP  IN  FIB  POLYMER: 

A  SEQUENTIAL  PROCESS?^^ 

8.1  OVERVIEW 

Oxygen  quenching  of  PtTFPP  in  FIB  polymer  is  analyzed  using  a  discrete  exponential  model  (DEM), 
a  Forster  decay  model  (FDM),  and  a  maximum  entropy  model  (MEM).  The  FDM,  which  has  been 
advocated  for  a  similar  system  to  counter  presumed  nonphysical  dependence  of  the  DEM  amplitudes  on 
quencher  density,  is  discounted  for  the  present  system.  Instead,  the  MEM  is  found  to  support  the  DEM  in 
suggesting  that  the  decay-rate  distribution  is  bimodal,  with  fractional  contributions  that  depend  on 
quencher  density.  A  sequential  decay  mechanism  involving  a  to-be-identified  intermediate  state  is 
proposed  to  explain  this  behavior. 

8.2  INTRODUCTION 

The  study  of  the  fluorescence  lifetime  behavior  of  molecules  sequestered  in  a  host  matrix  continues  to 
be  the  subject  of  active  research.^’^^’*^’**’*^’**'^'"’  Since  about  1985,  it  has  been  clear  that  it  can  be  difficult 
to  ascertain,  by  analysis  of  the  fluorescence  decay  alone,  whether  the  underlying  lifetime  spectrum  is 
continuous  or  discrete  in  nature. Additional  insight  into  the  dynamics  may  be  obtained  by  studying  the 
quenching  behavior  of  such  systems.  88,89,105-110  particular,  the  fractional  contributions  of  discrete 
components  in  the  spectrum  should  be  independent  of  quencher  density  if  these  contributions  are  to 
reflect  a  fixed  partitioning  of  the  fluorophore  molecules  with  regard  to  quenching.  ’  ’  In  1995,  such 
considerations  led  Draxler  et  al.  to  advocate  a  Forster  decay  model  (FDM)  over  a  discrete  exponential 
model  (DEM)  for  oxygen  quenching  of  Ru(dpp)  in  a  number  of  polymers.**’*^  Here,  the  opposite  finding 
is  reported  for  oxygen  quenching  of  PtTFPP  in  FIB  polymer  [Pt(II)  meso-tetra  (Pentafluorophenyl) 
porphine  in  fluoro-isopropyl-butyl],  a  sprayable  paint  that  has  become  popular  for  measuring  air-pressure 
distributions  on  test  articles  in  wind  tunnels  by  optical  means. An  alternative  explanation, 
involving  sequential  decay,  is  proposed  for  the  observed  lifetime  behavior  of  this  system.  Results  of  the 
analysis  are  supported  by  a  maximum  entropy  model  (MEM).^^’^^’^®’^*’^°° 

8.3  EXPERIMENT 

The  paint  PtTFPP  in  FIB  was  developed  at  the  University  of  Washington. and  is  available  from 
Innovative  Scientific  Solutions,  Inc.  Painted  aluminum  coupons  were  mounted  in  a  calibration  chamber  in 
which  air  pressure,  low  humidity,  and  sample  temperature  were  precisely  controlled.  Samples  were 
illuminated  by  an  array  of  blue  EEDs  with  peak  emission  at  465  nm,  resulting  in  excitation  of  the  PtTFPP 
first  excited  singlet  state.  The  illumination  level  was  held  constant  to  better  than  1  percent.  Fluorescence 
from  the  triplet  state,  filtered  at  650  nm  using  a  40-nm  FWHM  bandpass  filter,  was  measured  by  a  fast 
photomultiplier  tube  (PMT).  The  EED  source  was  pulsed  at  4.95  kHz.  Nominal  pulse  width  was  10  ps, 
with  a  rise  time  of  approximately  0.5  ps.  Output  from  the  PMT  was  sampled  at  0.1 -ps  intervals,  from 
-0.5  ps  to  about  -l-lOO  ps  relative  to  the  trigger  input,  and  this  process  resulted  in  1002  samples  per  curve. 
Curves  were  summed  over  100  pulses  using  a  digital  oscilloscope.  Preprocessing  of  the  data  was 
performed  to  subtract  a  baseline  value  and  to  determine  the  precise  start  of  the  optical  pulse;  the  result 
was  fluorescence  response  curves,  S(t),  such  as  those  shown  in  Fig.  50.  Measurements  were  carried  out 
over  a  range  of  air  pressures  (from  1  to  100  kPa)  and  temperatures  (from  4  to  49°C). 


§§Prgviously  published  as  “Oxygen  Quenching  of  PtTFPP  in  FIB  Polymer:  a  Sequential  Process?”  in  Chemical 
Physics  Letters,  Vol.  394,  No.  1-3, 2004,  pp.  101-104. 
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8.4  ANALYSIS 

The  measured  fluorescence  response  curves,  S(t),  were  fit  to  delta-function  response  curves,  Ss(t), 
convolved  with  the  measured  pulse  shape;  also  taken  into  account  were  small  amounts  of  residual 
fluorescence  (<  4  percent)  from  prior  pulses.*'  For  the  DEM,  FDM,  and  MEM,  Ss(t)  is  presumed  to  be 
given,  respectively,  by  the  following  relations  (see  also  Section  7.3): 

2 

k=l 

Ss,FDM(0  =  Soexp(-at-bt''^),  (87) 

and 

100 

=  (88) 

fi=i 

The  DEM  and  FDM  have  four  and  three  fit  parameters,  respectively,  while  the  MEM  has  100.  (The 
MEM  decay  rates,  y^,  were  fixed  in  equal  increments,  Ay^,  along  a  log  axis  that  stretches  from  y  =  0.001 
ps''  to  y  =  10  ps’'.)  Smoothness  of  the  MEM  distribution  was  enforced  by  a  Shannon-Jaynes  entropy 
constraint,  following  Refs.  93,  98,  and  100.  All  analysis  was  based  on  least-squares  fitting  of  the 
measured  fluorescence  response  curves,  including  the  excitation  phase.  Further  details  of  the  experimental 
and  numerical  techniques  employed  may  be  found  in  Ref.  8 1 . 

8.5  RESULTS 

Figure  5 1  shows  selected  fit  residuals  for  the  three  models.  The  three-parameter  FDM  consistently 
produced  fits  that  were  equally  good  or  somewhat  better  than  the  four-parameter  DEM  fits,  while  the 
MEM  produced  the  best  fits  overall.  Relative  RMS  fit  errors,  s  (expressed  as  a  percentage  of  the  peak 
signal),  were  in  the  range  from  0.35  to  0.65  percent  for  all  data. 

Figure  53  shows  the  implied  decay-rate  distributions  for  the  data  from  Fig.  51.  For  the  FDM,  this 
distribution  is  obtained  by  taking  the  inverse  Eaplace  transform  of  Eq.  (87),  yielding^^ 


where  G(x)  =  x*  exp(-x^)  and  ny(y)  =  0  for  y  <  a.  For  the  MEM,  the  distribution  in  Fig.  53  is  given  by  the 
factors  in  Eq.  (88),  divided  by  the  widths,  Ay^,  of  the  y^  bins.^^  While  the  modality  of  the  DEM  and  the 
FDM  is  inherent  in  the  definition  of  the  presumed  delta-function  response  (bimodal  and  unimodal, 
respectively),  the  modality  of  the  MEM  is  not  determined  a  priori.  For  this  reason,  the  bimodal  nature  of 
the  MEM  distribution  may  be  taken  as  a  strong  indication  that  the  fluorescence  decay  is  governed  by  two 
time  constants,  or  at  least  by  two  distinct  distributions.  (As  established  previously,  the  widths  of  peaks  in 
a  MEM  spectrum  are  attributable,  in  part,  to  noise  on  the  data.^*)  Further  indication  that  the  FDM  is  not 
the  correct  model  for  the  present  system  lies  in  the  fact  that  the  fit  parameter  a  from  Eq.  (87)  is  found  to 
decrease  with  pressure  (a  roughly  linear  increase  would  be  expected’"),  while  the  fitted  value  b  from  Eq.. 
(87)  is  found  to  be  strongly  dependent  on  pressure,  instead  of  being  approximately  constant.*' 

On  the  other  hand,  while  the  results  of  the  DEM  produce  the  largest  fit  residuals,  they  are  remarkably 
consistent  with  results  from  the  MEM,  both  with  regard  to  the  positions  of  the  two  peaks  in  the 
distribution  (see  Fig.  56)  and  with  regard  to  the  fractional  contributions  of  the  two  peaks  to  the  total  signal 
(see  Fig.  57).  (For  the  MEM,  the  peak  positions  are  those  of  the  maxima  in  the  distributions,  while  the 
fractional  contributions  were  determined  by  dividing  the  MEM  distributions  at  the  minima  between  the 
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peaks;  the  data  in  Fig.  57  were  also  averaged  over  the  range  from  4  to  49°C,  because  the  fractional 
contributions  were  found  to  be  insensitive  to  temperature.) 

Two  features  stand  out  in  Figs.  56  and  57:  the  remarkable  linearity  of  the  Stem-Volmer  (SV)  behavior 
in  Fig.  56  and  the  strong  dependence  of  the  fluorescence  fractions  on  pressure  in  Fig.  57.  This  latter 
dependence  is  inconsistent  with  a  fixed  partitioning  of  the  fluorophore  molecules  with  respect  to 
quenching  and  is  precisely  the  type  of  behavior  that  led  Draxler  et  al.  to  advocate  the  FDM  over  the  DEM 
for  their  system.** 

8.6  PROPOSED  MODEL 

It  is  now  argued  that  the  behavior  shown  in  Fig.  57  is  consistent  with  a  sequential  decay  model,  in 
which  two  microsecond-type  processes  play  a  role.  This  concept  is  illustrated  in  Fig.  59,  in  which  two 
excited  fluorophore  states,  “1”  and  “2,”  are  populated  at  rates  pi  and  pi  via  the  short-lived  singlet  state, 
are  depopulated  at  rates  yi  and  72,  and  have  a  one-way  crossover  rate,  712.  Fluorescence  dynamics  of  this 
kind  are  found  in  donor-acceptor  systems,  in  which  excited-state  energy  is  transferred  from  one  type  of 
fluorophore  to  another, yielding  fluorescence  lifetimes  ti  =  1/(71  -l-  712)  and  X2  =  I/72.  Assume  now  that 
the  decay  rate  71  is  negligible  (i.e.,  71  =  0).  In  this  case  it  is  not  necessary  to  specify  the  relative 
efficiencies  with  which  fluorescence  from  states  1  and  2  is  detected,  and  the  fractional  contribution  of 
lifetime  ti  to  the  fluorescence  signal  is  given  by  an  expression  that  is  similar  to  that  which  follows  from 
Eqs.  (17)  and  (18)  in  Ref.  112,  with  p\  and  p2  playing  the  roles  of  the  initial  donor  and  acceptor 
populations,  71  and  72  playing  the  roles  of  the  donor  and  acceptor  decay  rates,  and  712  playing  the  role  of 
the  donor-acceptor  transfer  rate.  Thus: 

(90) 

Pi  1y2  -Yi2, 


First  Singlet  State  (SI) 


Figure  59.  Schematic  of  Proposed  Sequential  Decay  Process 

Now  assume  that  the  rates  712  and  72  are  given  by  the  solid  lines  in  Fig.  60a,  which  are  linear  fits  to  the 
MEM  results  in  Fig.  56.  Also  assume  that  p]  =  pi.  The  fluorescence  fraction, /i,  from  Eq.  (90)  and  its 
complement,  fi  =  \  -  f\  ,  are  then  given  by  the  solid  lines  in  Fig.  60b.  The  agreement  between  this 
predicted  behavior  and  the  measured  fractions  from  Fig.  57  is  remarkable  and  thus  argues  for  the  validity 
of  the  proposed  model. 
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a.  Peak  Positions 


P,  kPa 

b.  Fluorescence  Fractions 


Figure  60.  Peak  Positions  and  Contribution  Fractions  for  Sequential  Decay  Model 


8.7  DISCUSSION 

The  assumptions  yi  =  0  and  p\  =  p2  are  perhaps  overly  simplistic.  More  generally,  one  might  allow 
that  yi  ^  0,  that  yi  and  yi2  have  separate  SV  dependencies  (only  the  sum  dependence  is  measured  in  Fig. 
4a),  that  p]  ^  p2,  and  that  both  states  1  and  2  contribute  to  the  observed  fluorescence  signal,  with  relative 
efficiencies  r|i  and  ri2.  This  expanded  assumption  introduces  four  additional  parameters  into  the  model: 
one  unknown  SV  intercept,  one  unknown  SV  slope,  the  ratio  p\lp2,  and  the  ratio  r|i/r|2.  This  allows  for 
additional  combinations  of  parameters  that  can  produce  the  fluorescence  fractions  in  Fig.  60b.  Thus,  the 
fact  that  the  fractions  in  Fig.  60b  can  be  reproduced  without  the  additional  four  parameters  may  be  a 
coincidence.  Even  so,  the  proposed  model  makes  a  strong  case  for  the  presence  of  a  heretofore  unknown 
intermediate  state.  To  identify  the  nature  of  this  state,  more  research  is  required,  involving,  perhaps, 
spectrally  resolved  measurements  in  addition  to  temporally  resolved  ones.^**®  Among  the  possible 
implications  of  this  work  is  the  extension  of  luminescent  paint  techniques  from  measurement  of  pressure 
only  to  measurement  of  pressure  and  temperature  in  combination,  as  is  discussed  in  detail  in  Section  9.0. 
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9.0  OPTIMIZATION  OF  THREE-GATE  LIFETIME  PRESSURE-  AND 
TEMPERATURE-SENSITIVE  PAINT  MEASUREMENTS*** 


9.1  OVERVIEW 

A  theoretical  model  is  described  to  optimize  the  performance  of  a  lifetime-based  pressure-sensitive 
paint  system,  in  which  signals  are  collected  at  three  gates,  with  the  goal  of  performing  a  temperature- 
compensated  pressure  measurement,  or  a  measurement  in  which  both  pressure  and  temperature  are 
obtained.  While  the  three-gate  approach  eliminates  the  temperature  bias  of  a  conventional  two-gate 
measurement,  the  shot-noise-induced  error  in  a  three-gate  measurement  is  not  necessarily  better  than  the 
temperature  bias  in  a  two-gate  measurement;  in  fact,  it  is  potentially  much  worse.  This  effect  is  quantified 
by  considering  the  relative  temperature  sensitivities  of  the  paint  at  the  selected  gate  ratios,  both  for  results 
that  have  been  reported  to  date  in  the  literature  and  for  new  paints  that  might  be  developed  by  combining 
pressure-sensitive  and  temperature-sensitive  probes  in  a  single  binder. 

9.2  INTRODUCTION 

Pressure-sensitive  paint  has  established  itself  as  an  important  technique  for  measuring  pressure 
distribution  on  test  articles  in  aerodynamic  test  facilities,  especially  transonic  wind  tunnels.  The  principle 
of  the  technique  has  been  reviewed  extensively  elsewhere.'’^  In  its  most  common  form,  a  paint  containing 
luminescent  probe  molecules  is  applied  to  a  test  article,  which  then  is  illuminated  by  short-wavelength 
light  to  produce  a  pressure-dependent  emission  signal  at  a  longer  wavelength.  The  pressure  dependence 
typically  comes  about  because  of  radiationless  quenching  of  excited  probe  molecules  by  oxygen. 
Measurements  can  be  performed  using  continuous  or  pulsed  illumination.  In  the  latter  case,  images  of  the 
test  article  are  captured  using  digital  cameras  that  are  gated  with  respect  to  the  train  of  illumination 
pulses.  By  suitable  selection  of  the  pulse  repetition  frequency,  pulse  duration,  and  start  and  end  times  of 
the  camera  gates,  the  ratio  of  the  signals  collected  during  two  gates  becomes  a  function  of  pressure, 
thereby  forming  the  basis  for  a  pressure  measurement. 

Among  the  complications  of  such  PSP  measurements  is  the  fact  that  signal  ratios  usually  depend  on 
temperature  to  a  mild  degree. Variations  in  surface  temperature  thus  tend  to  produce  pressure 
measurement  errors.  Recent  work  on  lifetime-based  PSP  suggests  that,  given  a  suitable  paint, 
measurement  at  a  third  gate  might  be  added  to  compensate  for  this  unwanted  temperature  dependence  and 
to  perform  a  temperature  measurement  in  its  own  right.  In  Section  9.3  this  previous  work  is  reviewed  and 
cast  in  a  form  that  lends  itself  to  a  new  theoretical  approach.  This  new  approach  is  developed  in  Section 
9.4  to  guide  further  work  in  this  area.  In  Section  9.5  the  model  from  Section  9.4  is  applied  to  published 
data  on  current  paints,  with  the  goal  of  optimizing  three-gate  PSP  measurements  and  defining  critical 
parameters  for  a  possible  new  paint.  Section  9.6  concludes  with  a  discussion. 

9.3  REVIEW  OF  PREVIOUS  WORK 

Figure  61  shows  an  example  of  a  fluorescence  response,  S(t;P,T),  that  might  be  the  result  of  excitation 
of  the  paint  by  a  square  pulse.  The  exact  shape  of  the  response  (especially  that  of  the  decay  portion) 
depends  on  the  pressure,  P,  and  temperature,  T,  of  the  surface.*®’*'  High  pressures  produce  faster  decays 
than  do  lower  pressures.  To  a  lesser  extent,  high  temperatures  also  produce  faster  decays  than  do  lower 
temperatures.  Because  PSP  is  primarily  an  image-based  technique,  it  is  not  practical  to  measure  the 
detailed  shape  of  the  response  function,  except  in  calibration  measurements.  Instead,  lifetime-based  PSP 
measurements  rely  on  integrating  the  signal  under  limited  portions  of  the  response  function.  The  single¬ 
pulse  signal.  Si,  for  a  gate  with  index  i  may  thus  be  written  as 


***Previously  published  as  “Optimization  of  Three-Gate  Lifetime  Pressure-  and  Temperature-Sensitive  Paint 
Measurements,”  AIAA  Paper  2004-2190,  24th  AIAA  Aerodynamic  Measurement  Technology  and  Ground  Testing 
Conference,  Portland,  OR,  June  28- July  1,  2004. 
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s,  =  V'S{t-P,T)&t,  (91) 

Jtl2 

where  tn  and  ta  denote  the  start  and  end  times  of  gate  i  (see  Fig.  61).  If  Ni  is  understood  to  represent  the 
number  of  pulses  over  which  the  gate-/  signal  is  integrated  and  N  to  represent  the  total  number  of  pulses 
(i.e.,  N  =  N\+  N2  for  a  two-gate  measurement  and  N  =  Ni+  N2+  N},  for  a  three-gate  measurement),  the 
pulse-integrated  signal  for  gate  /  may  then  be  written  as 

Si^NfiSi,  (92) 

where =  NilN.  That  is,  ft  represents  the  fraction  of  the  total  measurement  that  is  devoted  to  collecting  the 
signal  at  the  /-th  gate.  (As  shown  by  Bell,  it  should  not  be  considered  a  foregone  conclusion  that  the  same 
number  of  pulses  should  be  devoted  to  each  gate.^) 


Fluorescence 

Response 


Figure  61.  Schematic  of  Gated  Measurement 

Based  on  Eq.  (92),  the  ratio  between  two  signals,  5,  and  5),  acquired  at  separate  gates  /  and  j,  may  be 
written  as 


Rh 


Sj  _  fjSj  _fi 
Sj  ■  fjSj  -  fj  ’ 


(93) 


where  ry  may  be  thought  of  as  the  intrinsic  single-pulse  ratio  for  the  two  gates.  Figure  62  shows  an 
idealized  dependence  of  this  ratio  on  pressure  and  temperature,  with  5,  measured  during  or  immediately 
following  the  excitation  pulse,  and  5,  measured  after  some  delay. 


Figure  62,  Dependence  of  Signal  Ratio  on  Pressure  and  Temperature 


If  the  paint  response  can  be  characterized  by  a  single  time  constant,  x,  the  ratio  Ry  may  be  written 
explicitly  in  terms  of  this  time  constant.  Given  these  conditions,  it  is  possible  to  define  the  relative 
sensitivity  {dx/x)/(dRy/Rij).  Several  authors  have  followed  this  approach  to  optimize  the  operational 
parameters  of  a  two-gate  PSP  measurement. Flere  a  more  general  approach  is  employed  that  allows 
for  the  presence  of  multiple  time  constants  in  the  fluorescence  response.  In  this  generalized  approach,  the 
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intrinsic  pulse  ratio  from  Eq.  (93)  at  a  pressure  Pq  +  AP  and  a  temperature  Tq  +  AT  near  an  operating 
point  {Po,To)  is  written  as 

nj=f^[  +  Pij^  +  Hj^T],  (94) 

where  is  the  ratio  at  the  operating  point  and  pij  and  ty  are  the  relative  partial  derivatives  of  this  ratio 
with  respect  to  pressure  and  temperature,  also  at  the  operating  point. 

Equation  (94)  is  not  meant  to  describe  accurately  the  calibration  equation  for  a  PSP  measurement 
over  the  entire  range  of  pressures  and  temperatures  of  interest.  This  generally  requires  the  inclusion  of 
quadratic  and  cubic  terms  in  addition  to  the  linear  terms  from  Eq.  (94).  Rather,  the  linear  form  of  Eq.  (94) 
is  used  in  Sections  9.4  and  9.5  to  approach  the  optimization  of  a  two-gate  or  three-gate  PSP  measurement 
near  an  operating  point  (EcTo). 

Table  1 1  lists  values  of  py,  and  ty  that  were  calculated  from  published  data  from  authors  who 
have  reported  early  results  on  combined  pressure  and  temperature  measurements  using  a  lifetime-based 
PSP  system.  These  are  Mitsuo  et  al.,^^  Watkins  et  al.,’^  Davies  et  al.,"^  and  Hradil  et  al.’°  Details  on 
paints  and  gating  parameters  are  noted  in  the  table  and  are  discussed  further  in  Section  9.5.  Of  the  authors 
listed  in  Table  11,  only  Hradil  et  al.  deliberately  introduced  a  second  time  constant  into  the  paint,  by 
adding  a  second  probe  molecule  (a  long-lived  phosphor  with  a  3-ms  time  constant). 

Table  11.  Paint  Parameters  at  1  atm  and  25°C 


Authors 

Paint 

Gate  1 

Gate  2/3 

r..(o) 

'u 

Pij,  atm  ' 

tii.^C-' 

i1ii,%/‘’C 

Mitsuo  at  al. 

PtTFPP  in  FEM 

0  to  0.8  ps 

1 2  to  1 9  ps 

0.70 

1.29 

0.0370 

-2.90 

0  to  0.8  ps 

30  to  83  |ts 

0.70 

1.55 

0.0670 

-4.30 

Watkins  et  al. 

PtTFPP  in  FEM 

0  to  1  ps 

12  to  20  ps 

2.10 

0.99 

0.0110 

-1.14 

0  to  1  ps 

30  to  83  |ts 

3.00 

1.23 

0.0150 

-1.20 

Watkins  et  al. 

PtTFPP  in  FIB 

0  to  1  ps 

12  to  20  ps 

1.76 

0.74 

0.0033 

-0.45 

0  to  1  ps 

30  to  83  |ts 

2.34 

0.99 

0.0048 

-0.49 

Davies  et  al. 

Porphyrin  in  Siloxane 

0  to  0.5t 

t:  to  2x 

1.69 

0.70 

0.0037 

-0.53 

Hradil  et  al. 

Ru{dpp)/MFG  in  MTEOS 

0  to  0.65  ps 

1 .3  to  2.6  (ts 

0.73 

0.08 

0.0022 

-2.61 

0  to  0.65  ITS 

3.5  to  7.1  ms 

0.014 

-0.05 

0.0002 

3.65 

From  Eq.  (94)  it  follows  that  the  temperature  sensitivity  of  a  PSP  measurement  can  be  expressed  as 


AT  Pij  ' 


(95) 


This  sensitivity  is  typically  quoted  as  percent  pressure  per  unit  temperature  at  some  operating 
pressure  Po  and  is  expressed  here  as 


1  AP  1  Hi 

rtij  = - X 100%  = - ^  X 100%  . 

'y  Pq  AT  Po  Pij 


(96) 


Values  of  this  relative  temperature  sensitivity,  Py,  are  included  in  Table  1 1  and  are  typically  on  the 
order  of  one  percent  per  degree  Celsius  at  one  atmosphere.  Such  values  are  also  typical  for  intensity- 
based  PSP  measurements,  in  which  signal  ratios  are  calculated  between  reference  (“wind-off’)  and  run 
(“wind-on”)  conditions  under  constant  illumination. ' 


9,4  THREE-GATE  VS,  TWO-GATE  MEASUREMENTS 

The  goal  in  this  section  is  to  take  a  conceptual  look  at  how  a  temperature  error  produces  a  pressure 
error  in  a  two-gate  PSP  measurement  and  how  this  error  might  be  reduced  by  a  three-gate  measurement. 
Eet  ASi  and  A5,  denote  errors  in  the  signals  Si  and  Sj  attributable  to  noise  or  systematic  effects  and  let  AP 
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and  Ar  denote  the  corresponding  errors  in  pressure  and  temperature.  From  Eqs.  (93)  and  (94),  it  follows 
that  these  errors  are  related  by 


Si 


AS/ 


(97) 


From  this  it  follows  [cf.  Eq.  (95)]  that  a  two-gate  PSP  measurement  has  a  bias  error  given  by 

P\2 


(98) 


where  AT  is  the  (generally  unknown)  amount  by  which  the  assumed  temperature  differs  from  the  actual 
temperature.  The  corresponding  noise-induced  error  is  calculated  by  assuming  (cf  Belf  and  Goss^*)  that 
the  signals  A,  and  Sj  are  shot-noise  limited  (the  best  possible  scenario),  with  a  relative  variance  given  by 

{{ASi! =  (j  likewise).  (99) 

\  /  kSi  kNfiSi 


FI  ere  k  is  the  number  of  detected  photons  per  unit  signal.  As  shown  in  Section  9.7,  the  root-mean- 
square  (RMS)  pressure  error  in  a  two-gate  measurement  may  then  be  written  as 


(2)  _ 
P  ~ 


\P12 


4kN 


(100) 


where  N  is  the  total  number  of  pulses  among  the  two  gates,  j3i2  is  the  relative  partial  pressure  sensitivity  at 
the  operating  point,  and  si  and  S2  are  the  single -pulse  signals  from  Eq.  (91).  The  result  from  Eq.  (100) 
assumes  that  the  pulse  allocation  fractions  f\  andy^  have  been  selected  optimally. 

In  a  three-gate  PSP  measurement,  the  bias  error  attributable  to  unknown  temperature  is  removed  (i.e., 
AP^^^  =  0),  and  the  RMS  pressure  error  becomes  (see  again  Section  9.7): 


^  I  P12T3  “P13T2  I 

Here  p\2  and  tn  are  the  relative  partial  pressure  and  temperature  sensitivities  for  the  signal  ratio  rn, 
and  pii  and  tu  are  the  relative  partial  pressure  and  temperature  sensitivities  for  the  signal  ratio  r^.  As  is 
the  case  for  the  two-gate  result  from  Eq.(lOO),  the  result  from  Eq.  (101)  assumes  that  the  pulse  allocation 
fractions  (in  this  case,/i,y^,  and/^)  have  been  selected  optimally. 

Before  the  above  errors  are  evaluated  numerically,  it  is  instructive  to  review  graphically  how  the 
errors  arise  in  a  P-T  diagram.  To  this  effect,  the  heavy,  solid  line  in  Fig.  63  represents  a  contour  of  equal 
signal  ratios  P 12  in  a  two-gate  measurement  as  a  function  of  pressure  and  temperature  for  a  noise-free 
measurement.  The  actual  values  of  pressure  and  temperature  are  represented  by  the  operating  point 
(RcTo).  In  the  two-gate  case  from  Fig.  63,  the  result  of  a  temperature  error  AT  is  to  cause  a  pressure  bias 
error  AP^^'  that  depends  on  the  relative  temperature  sensitivity  pn  [see  Eqs.  (98)  and  (96)].  Because  of 
shot  noise,  the  equi-ratio  contour  is  broadened,  as  shown  by  the  gray  band  in  Fig.  63,  and  the  result  is  a 
corresponding  random  component  of  pressure,  with  a  magnitude  given  by  Eq.  (100). 

By  contrast,  there  is  no  bias  error  in  the  three-gate  case  depicted  in  Fig.  64,  by  virtue  of  the  fact  that 
the  operating  point  (P’0,  Po)  is  defined  by  the  intersection  of  two  equi-ratio  contours,  Ru  and  Ru. 
However,  a  potentially  much  larger  random  error,  results,  depending  on  the  shape  of  the  overlap 
region  of  the  associated  noise  bands  in  Fig.  64.  From  Fig.  64  it  is  clear  that  the  two  ratio  contours  in  a 
three-gate  measurement  must  intersect  at  an  appreciable  angle  in  order  not  to  stretch  the  diamond-shaped 


|t3  -T2I  ^  |t3|  ^  IT2I 
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(101) 
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Figure  63.  Error  Sources  in  Two-Gate  PSP 


overlap  region  significantly.  This  effect  is  controlled  by  the  leading  denominator  in  Eq.  (101),  which 
deviates  from  zero  precisely  to  the  degree  by  which  the  temperature  sensitivities  pn  and  pn  of  the  two- 
gate  ratios  differ  from  each  other. 

Whether  the  three-gate  RMS  error  from  Eq.  (101)  is  smaller  than  the  two-gate  bias  error  from  Eq. 
(98)  depends  on  several  factors — specifically,  the  temperature  sensitivity  of  the  two-gate  measurement 
(which  can  be  calibrated),  the  temperature  error  in  the  two-gate  scheme  (which  is  generally  unknown), 
and  the  number  of  pulses  over  which  the  three-gate  measurement  is  performed. 

Comparison  of  the  three-gate  RMS  error  from  Eq.  (101)  to  the  two-gate  RMS  error  from  Eq.  (100)  is 
possible  if  1)  the  same  number  of  pulses,  N,  is  used  in  both  cases;  2)  the  single-pulse  signals  si,  S2,  and  S3 
in  Eqs.  (100)  and  (101)  are  roughly  balanced  in  the  sense  that  (for  some  common  so)  si  =  52  =  so/2  in  Eq. 
(100)  and  si  =  52  =  53  =  5o/3  in  Eq.  (101);  and  3)  tu  >  tn  >  0.  In  terms  of  the  temperature  sensitivities  pi2 
and  p  13  for  the  two  ratios  from  Eq.  (96),  it  then  follows  that 


713 

712  -713 


(102) 


This  result  suggests  again  that,  barring  large  temperature  errors  in  the  two-gate  scheme,  a  necessary 
condition  for  a  three-gate  measurement  to  outperform  a  two-gate  measurement  is  that  the  relative 
temperature  sensitivities  pn  and  pn  in  the  three-gate  scheme  must  be  substantially  different  from  each 
other. 


9,5  PAINT  OPTIMIZATION 

The  requirement  that  the  relative  temperature  sensitivities  pi2  and  p  13  in  a  three-gate  measurement 
scheme  be  substantially  different  from  each  other  appears  to  be  marginally  the  case  for  the  work  on 
PtTFPP  in  FEM  by  Mitsuo  et  al.^^  in  Table  11,  with  p  12  — 2.9  percent  per  degree  Celsius  and  pi3  — 4.3 
percent  per  degree  Celsius  at  ambient  conditions.  According  to  Eq.  (102),  this  gives  a  three-gate  RMS 
error  that  is  nominally  three  times  as  large  as  the  RMS  error  in  a  two-gate  measurement. 

Quite  different  results  are  obtained  by  Watkins  et  al.^^  for  a  similar  system  (PtTFPP  in  FEM  also)  and 
equal  gate  times,  i.e.,  p  12  =  -1.14  percent  per  degree  Celsius  and  po  =  -1.20  percent  per  degree  Celsius 
in  Table  11.  The  fact  that  these  relative  sensitivities  are  equal  within  the  accuracy  of  the  data  suggests 
that  a  useful  three-gate  measurement  is  not  possible. 

Eikewise,  for  PtTFPP  in  FIB,  Watkins  et  al.^^  report  results  that  yield  relative  sensitivities  of  pi2  = 
-0.45  percent  per  degree  Celsius  and  p  13  =  -0.49  percent  per  degree  Celsius,  as  shown  in  Table  11.  These 
values  also  seem  to  preclude  a  useful  three-gate  measurement,  despite  efforts  by  the  authors  of  Ref  1 3  to 
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use  the  paint  in  this  manner.  (Unpublished  analysis  by  the  present  author  on  PtTFPP  on  FIB  reaches  the 
same  conclusion:  that  PtTFPP  in  FIB  does  not  lend  itself  to  a  three-gate  measurement  scheme.) 

It  is  not  clear  why  Mitsuo  et  al.  and  Watkins  et  al.  obtain  different  results  for  PtTFPP  in  FEM,  but  the 
results  for  PtTFPP  in  FIB  from  Ref  1 3  are  not  too  surprising,  in  light  of  the  fact  that  this  paint  has  been 
described  as  “ideal”  by  its  developers,  in  the  sense  that  its  temperature  sensitivity  is  independent  of 
pressure.^  Moreover,  PtTFPP  in  FIB  was  developed  with  minimum  temperature  sensitivity  in  mind. 

The  authors  of  Ref  113  have  also  claimed  success  with  a  multigate  scheme  for  measuring  pressure 
and  temperature,  but  they  have  not  formally  published  their  results  beyond  the  early  report  from  Ref  113. 

Regarding  results  published  to  date,  this  leaves  the  work  of  Hradil  et  al.  from  Ref  1 0,  in  which  a 
temperature  probe  (the  long-lived  phosphor  MFG)  was  added  to  a  conventional  PSP  paint,  Ru(dpp)  in 
MTEOS.  Although  Hradil  et  al.  do  not  consider  a  gating  scheme  per  se,  such  a  system  may  be  analyzed 
by  assuming  that  the  fluorescence  response  can  be  characterized  by  two  decay  rates,  yi  and  72,  with 
relative  contributions  in  the  proportion  x  :  (1-x).  The  response  of  such  a  system  to  delta- function 
excitation  may  be  written  as 

S(t)  =  xe~^‘^  +(l-x)e~^2^  .  (103) 

With  Eq.  (91)  it  follows  that  the  resulting  single -pulse  signal  for  a  gate  i  is  then  given  by 

^  j  +  (Izf)  j ,  (104) 

The  pressure  and  temperature  sensitivities  of  the  signal  ratio  r,j  =  s/sj  may  now  be  written  in  terms  of 
the  partial  derivates  of  the  rates  yi  and  72  with  respect  to  F  and  T  as 


and 


1  Ssj 

1  dsi 

1  dSj^ 

dn 

ry  dP 

Sj  dn^ 

dP 

^72 

sj  dr2 ) 

dP 

1  dsi 

1  dsA 
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1 

dsi 

1  dSj'] 

3/2 

^  dri 

sj  dnj 

dT 

8/2 

sj  drzj 

dT 

(105) 


(106) 


Given  values  of  yi,  dr^\ldP,  and  dr^\ldT  for  the  first  component  and  values  72,  d^iIdP,  and  for  the 

second  component,  the  optimization  of  the  three-gate  measurement  error  from  Eq.  (101)  may  now  be  cast 
as  the  problem  of  establishing  the  optimum  ratio  x  in  Eq.  (103),  along  with  finding  optimum  gate  times 
t\\,  t\2,  hi,  tn,  hu  h2,  optimum  pulse  allocation  fractions  f\,f2,f2,  an  optimal  pulse  duration,  tp;  and  an 
optimal  pulse  repetition  period,  4-  A  representative  (nonoptimized)  result  is  quoted  in  Table  1 1  for  Hradil 
et  al.,  assuming  that  x  =  0.5,  yielding  =  -2.61  percent  per  degree  Celsius  and  po  =  3.65  percent  per 
degree  Celsius.  These  values  produce  a  factor  0.6  when  substituted  into  Eq.  (102),  suggesting  that  the 
implied  three-gate  temperature  measurement  not  only  is  free  from  a  temperature-bias  error  but  also  has  an 
RMS  error  that  is  on  the  same  order  as  the  RMS  error  in  a  two-gate  measurement. 

Incidentally,  it  follows  from  Eqs.  (103)  through  (106)  that  fluorescence  decay  that  can  be 
characterized  by  a  single  time  constant,  t,  does  not  lend  itself  to  a  three-gate  scheme.  This  can  be  seen  by 
setting  X  =  1  in  Eq.  (103),  in  which  case  the  ratio  tijlpij  may  be  reduced  to 
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h^_d^  jdr 

Pij  “sr/ ap' 

This  result  is  independent  of  the  indices  i  and  j  and  thus  of  the  gating  scheme. 


(107) 


9.6  DISCUSSION 

The  results  from  Section  9.5  demonstrate  that  a  pressure-sensitive  paint  that  is  specifically  designed 
to  be  sensitive  to  temperature  (namely,  by  deliberate  addition  of  a  temperature  probe)  may  lend  itself 
better  to  a  three-gate  measurement  scheme  than  does  a  paint  in  which  only  residual  temperature 
dependence  is  present  (i.e.,  unintentional  temperature  dependence  that  could  not  be  removed  during  the 
design  of  the  paint).  In  particular,  the  theoretical  model  that  is  developed  here  indicates  that  a  critical 
performance  parameter  for  a  three-gate  measurement  scheme  is  the  degree  to  which  the  relative 
temperature  sensitivities  of  the  two  signal  ratios  differ  from  each  other. 

Other  considerations  have  yet  to  be  considered  in  full.  For  example,  the  temperature  probe  advocated 
by  Hradil  et  al.  in  Ref  13  has  a  lifetime  that  is  several  orders  of  magnitude  greater  than  that  of  the 
pressure  probe.  This  means  that  the  repetition  rate  of  the  excitation  source  has  to  be  reduced  significantly, 
and  that,  as  a  result,  measurement  time  will  be  increased.  A  temperature  probe  with  a  significantly  shorter 
lifetime  would  alleviate  this  concern. 

Other  practical  considerations  will  also  have  to  be  considered,  such  as  the  spatial  homogeneity  of  the 
paint,  particularly  in  the  case  that  different  probe  molecules  are  mixed  together  in  a  single  binder. 

Finally,  optimization  of  a  PSP  measurement  system  has  to  be  considered  not  only  for  a  fixed 
operating  pressure  and  temperature,  but  also  for  the  entire  range  of  pressures  and  temperatures  that  are  of 
interest  in  the  measurement. 

In  summary,  considerable  further  work  is  required  to  optimize  three-gate  PSP  measurements.  The 
theoretical  model  developed  in  this  paper  should  serve  as  a  useful  tool  for  such  further  work. 


9.7  APPENDIX  9A:  DERIVATION  OF  TWO-GATE  AND  THREE-GATE  SHOT-NOISE 
ERRORS 

When  the  temperature-induced  bias  error  for  a  two-gate  PSP  measurement  is  ignored,  it  follows  from 
Eq.  (107)  that  the  pressure  error  is  given  by 


PU 
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From  this  and  Eq.  (99)  it  follows  that  the  resulting  RMS  pressure  error  is  given  by 
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With  regard  to  the  pulse  allocation  fractions  f\  andy),  this  RMS  error  is  minimized  by  the  choice 


v^l  v^2 


(110) 


That  is,  most  of  the  measurement  time  should  be  devoted  to  the  gate  with  the  smallest  signal  per 
pulse.  From  Eqs.  (109)  and  (110),  the  result  from  Eq.  (100)  is  obtained. 


85 


AEDC-TR-06-6 


When  a  third  signal,  S3,  is  measured,  a  second  signal  ratio  may  be  defined,  say,  =  5'i/5'3  = 
with  associated  relative  partial  pressure  and  temperature  sensitivities  p\3  and  tn.  A  result  similar  to  Eq. 
(97)  may  be  written  down  for  this  second  ratio.  When  combined  with  Eq.  (97)  for  the  first  ratio,  Rn,  a 
linear  system  of  equations  in  terms  of  the  pressure  error,  AP,  and  the  temperature  error,  ST,  results.  This 
system  may  be  solved  to  yield  the  following  errors  for  pressure  and  temperature  in  terms  of  the  errors  in 
the  signals  5i,  S2,  and  S3. 
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It  can  be  shown  that  the  results  from  Eqs.  (Ill)  and  (112)  are  independent  of  the  ordering  of  the  gates 
in  the  ratios  (e.g.,  S2/S\  vs.  5'i/5'2)  and  independent  of  the  choice  of  common  gate  (e.g.,  S2  common  vs.  Si 
common,  as  is  assumed  here). 

When  the  signals  are  shot-noise  limited  according  to  Eq.  (99),  the  three-gate  RMS  pressure  error 
becomes 
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A  similar  result  can  be  derived  for  the  predicted  RMS  temperature  error.  In  general,  choices  of  gate 
times  and  pulse  allocation  fractions  that  optimize  the  pressure  error  from  Eq.  (113)  will  not  optimize  the 
corresponding  temperature  error,  and  vice  versa.  Here,  the  focus  is  on  the  pressure  measurement  error, 
which  is  minimized  with  respect  to  the  pulse  allocation  fractions  f\,f2,  and/)  when  [cf  Eq.  (110)]: 
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For  this  choice  of  pulse  allocation  fractions,  the  RMS  pressure  error  from  Eq.  (113)  may  be  written  in 
the  form  of  Eq.  (101). 
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